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ABSTRACT 


I  iO  * '  *-f 


The  objective  of  this  project  has  been  tO'design>an  Improved 
three-axis  stability  augmentation  system  (SAS)  for  the  AV-8B 
Advanced  Harrier  VSTOL  aircraft  using  microprocessor- based 

fry 

digital  control.  Th*  research^ focuses  on  improving  the  handling 
qualities  of  the  airplane  through  SAS  redesign  in  the  low  speed 
flight  regime.  Particular  attention  is  paid  to  the  so-called 
^weather -cocking1*  instability  encountered  in  transition 

SiO/'-p  r  O  rri) 

(hover^  conventional)  fl ight.  Until  it*  recently*' $here 

ij-torntc's 

has  beeii/a  dearth  of  information  about  the/f light  characteristics 
of  the  Harrier.  A  major  breakthrough  in  this  field  was  achieved 
by 'the  development  of  MCAIR's  X22A  AV-8B  mathematical  model, 
which  yielded  a  set  of  linearized  stability  derivatives  for  the 

A 

7 

aircraft.  The  first  step  toward  the  improvement  of  the  AV-8B  SAS 
requires  the  trt+Hzation^of  these  coefficients  in  the  development 

of  an  analog/hybrid  model  of  both  the  aircraft  itself  and  of  the 

f \X_sjls 

unmodified  SAS.  The  controller  design  employ^ digital  state  feed¬ 
back  control  to  relocate  the  system  closed  loop  poles.  The  coip  . 
elusion  reached  isHhat  this  method  of  control  represents  a  valid 
approach  to  the  final  solution  of  the  Harrier's  stability  and 

M  9  y  ‘ 

control  problems,  A  long  range  goal  of  this  research  is  that  this 

Cjvv i  ^  li  f— 

controller  design  concept,. be  applied  to  future  aerospace  vehicles. 


Accession  Fcr 


ITT  1 3  GiviAil 
DDC  TAB 
Unannounced 
Justification 


PREFACE 


_P i3trivut lor/ 


This  investigation  was  conducted  at  the  United  States 
Naval  Academy,  Annapolis,  Maryland  under  the  Trident  Scholar^ 
Program.  This  publication  is  the  final  technical  report. 


and  it  includes  results  through  30  May  1979.  An  oral  presenta¬ 
tion  was  made  on  24  April  1979  to  RADM  William  P.  Lawrence, 
Superintendent,  U.  S.  Naval  Academy  and  other  distinguished 
guests  at  the  annual  Trident  Scholar's  Dinner.  The  research 
was  conducted  in  the  Department  of  Weapons  and  Systems 
Engineering,  and  in  cooperation  with  the  Department  of 
Aerospace  Engineering. 
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special 


The  principal  investigator  for  the  study  was  myself. 
Midshipman  Robert  V.  Walters.  The  faculty  advisors  for  the 
project  were  Associate  Professor  E.  E.  Mitchell  (principal 
advisor)  and  Assistant  Professor  K.  A.  Knowles  of  the  Systems 
Engineering  Department,  and  CDR  M.  D.  Hewett,  Chairman  of 
the  Aerospace  Engineering  Department.  The  assistance, 
insistence,  and  perseverance  that  Professor  Mitchell,  Professor 
Knowles,  and  Commander  Hewett  have  shown  has  been  a  major 
catalyst  in  the  research  effort.  Considering  the  tremendous 
amount  of  engineering  expertise  and  experience  that  these 
men  brought  to  bear  on  the  problem,  it  is  not  surprising 
that  such  satisfying  results  were  achieved. 

Robert  V.  Walters 
Midshipman  First  Class, 
United  States  Navy 
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I  worry  that  the  existence  of  a  practical  VSTOL 
aircraft  today  (the  Harrier,  in  its  operational 
A  and  forthcoming  B  versions)  is  being  ignored 
in  favor  of  the  long-range  development  of  a 
"high-performance"  VSTOL  culminating  in  the 
1990s.  Even  this  development  has  been  slowed 
by  the  Defense  Department,  stating  that  before 
further  funds  are  devoted  to  VSTOL  it  must  be 
proven  superior  to  conventional  aircraft.  In 
my  judgment,  if  it  were  supported  at  this  time, 
an  effective  VSTOL  aircraft  could  be  developed 
as  early  as  1986.  The  Soviets  are  doing  it _ 


ADM  Elmo  R.  Zumwalt,  USN  (Ret) 
May  1979,  following  Department 
of  Defense  cancellation  of 
the  AV-?B  program. 

[17,  p.  104] 
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CHAPTER  1 


INTRODUCTION:  JUSTIFICATION  AND 
ORGANIZATION  OF  THE  REPORT 

One  of  the  major  controversies  facing  naval  policy 
makers  is  the  issue  of  the  suitability  and  effectiveness  of 
VSTOL  aircraft  in  a  ->ea-gow:,  combatant  environment.  Major 
policy  decisions  on  the  future  of  VSTOL  aircraft  are  being 
formulated  based  on  the  Marine  Corps  experience  with  the 
Hawker  Siddeley  AV-8A  Harrier.  The  success  or  failure  of 
the  AV-8B  now  in  development  will  further  influence  the 
Navy's  future  VSTOL  plans. 

The  AV-8A  has  experienced  a  very  high  accident  rate  in 
the  Marine  Corps.  In  the  first  six  and  one-half  years  of  its 
use  in  Marine  Corps  operation,  the  Harrier  logged  twenty- 
eight  crashes  —  nearly  or;  fifth  of  the  fleet  —  over  twenty 
of  which  resulted  in  the  total  destruction  of  the  aircraft. 

The  aircraft's  problems  mirror  the  present  state  of  the  art 
in  VSTOL  technology.  Present  technology  is  capable  of 
producing  engines  of  sufficient  thrust  at  low  weight  to 
provide  the  thrust-to-weight  ratios  demanded  in  VSTOL  airplane 
design.  Special  VSTOL  airframe  structural  and  aerodynamic 
considerations  are  also  well  within  present  technological 
capabilities.  Stability,  control,  and  handling  qualities 


problems,  however,  persist.  They  are  the  major  threat  to 
the  success  of  the  VSTOL  concept.  The  poor  transition-regime 
handling  qualities  of  the  AV-8A  due  to  inadequate  control 
authority  and  inadequately  designed  stability  augmentation 
are  well  documented  and  have  contributed  to  the  poor  reputation 
of  the  airplane  as  a  difficult,  unforgiving  airplane  to  fly. 

A  strong  research  effort  in  stability  augmentation  design 
is  required  to  insure  the  success  of  the  AV-8B  and  future 
VSTOl  efforts.  In  addition  to  the  motivation  resulting  from 
the  political  exigencies  surrounding  the  Harrier  program,  the 
strong  current  interest  in  the  study  of  digital  flight  control 
provides  further  justification  for  the  research. 

With  these  thoughts  in  mind,  it  was  decided  to  conduct 
an  investigation  of  the  flight  characteristics  of  the  AV-8B 
Advanced  Harrier,  with  the  intent  of  improving  the  effectiveness 
of  its  stability  augmentation  system  (SAS)  in  the  low  speed 
flight  regime.  The  goal  of  the  study  was  the  realization 
of  a  microprocessor-compatible  compensator  that  would  stabilize 
the  Harrier.  It  was  planned  that  the  new  system  would  be 
capable  of  functioning  at  two  levels  of  control  authority. 

The  "normal"  mode  would  be  designed  to  improve  the  overall 
handling  qualities  of  the  aircraft,  but  still  allow  the  pilot 
sufficient  control  to  meet  the  demanding  maneuverability 


requirements  of  a  modern  fighter/attack  aircraft.  The 
"recovery"  mode  would  essentially  take  control  of  the  aircraft, 
seeking,  in  spite  of  any  external  disturbances,  to  return 
the  aircraft  to  some  reference  "trim"  condition.  The  digital 
state  variable  control  1 2"  ..ha..  s.;.te  has  become  the  fruit  of 
the  study  lives  up  ^ar.y  u.  .  -tse  initial  expectations. 

This  report  documents  the  year-long  research  effort 
leading  up  to  the  final  control ier  design.  The  organization 
of  the  report  is  as  follows. 

The  body  of  the  report  is  divided  into  ten  chapters. 

The  first  three  chapters  are  composed  primarily  of  background 
information.  Chapter  2  describes  the  subject  of  the  investi¬ 
gation,  the  Harrier.  Particular  attention  is  paid  to  the 
existing  flight  control  and  stability  augmentation  systems. 

The  reasons  for  the  Harrier’s  unstable  transition  flight 
characteristics  are  discussed  in  Chapter  3.  The  next  two 
chapters  concern  the  formulation  of  a  mathematical  model 
for  the  Harrier.  Chapter  4  presents  a  derivation  of  the  VSTOL 
equations  of  motion,  while  Chapter  5  discusses  the  methodology 
for  the  experimental  determination  of  the  coefficients. 

Chapters  6  and  7  contrast  the  controller  design  of  the  present 
conventional  SAS  and  that  of  the  digital  state  feedback  device. 


A  description  of  the  digital  and  the  analog/hybrid  computer 
simulations  is  contained  in  Chapter  8.  Finally,  Chapter  9 
treats  the  results  of  the  investigation,  and  conclusions  are 
presented  in  Chapter  10.  Thirteen  Appendices  document  the 
computer  programming  necessary  in  the  research,  and  contain, 
in  the  form  of  computer-generated  time  histories,  the  results 
of  the  project. 


CHAPTER  2 
THE  HARRIER 

The  genesis  of  the  design  of  the  Harrier's  basic 
airframe  is  unique  in  +har  the-  entire  aircraft  was  quite 
literally  built  around  its  propulsion  plant  -  the  Rolls 
Royce  Pegasus  engine.  The  Pegasus  is  an  axial  flow  twin 
spool  turbo  fan  engine,  and  in  this  respect  is  similar  in 
principle  to  the  engines  powering  many  of  our  modern  airliners. 
The  requirements  of  VSTOL  flight  necessitate  several  departures 
from  the  norm,  however.  The  most  important  of  these  changes 
concerns  the  ducting  of  the  jet  exhaust.  Instead  of  the 
conventional  configuration  of  a  single  nozzle  located  in  line 
with  the  turbine  shafts,  the  Pegasus  uses  four  rotatable 
nozzles  located  symmetrica Hy  on  the  sides  of  the  engine 
(and  protruding  from  the  r>;<-laje  directly  under  the  shoulder- 
mounted  wings).  This  enqine-airframe configuration  gives  the 
Hawker  Siddeley  aircraft  unique  capabilities.  The  AV-8 
Harrier  has  the  ability  to  "vector'1  its  thrust,  i.e. ,  to  change 
both  the  direction  and  length  (magnitude)  of  the  thrust  vector. 
As  will  be  explained  in  a  moment,  thrust  vectoring  is  the 
Harrier's  means  of  achieving  VSTOL  capability.  Another 
unique  feature  is  that  the  two  compressor/turbine  sections 
of  the  Pegasus  engine  are  contrarotating.  By  spinning  the 
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blading  In  opposite  directions,  the  gyro-coupling  effects 
that  had  the  plagued  many  previous  VSTOL  designs  have  been 
all  but  eliminated.  An  additional  performance  criteria  met 
by  the  Pegasus -Harrier  package  results  from  the  location 
of  the  exhaust  nozzles  such  that  the  net  thrust  vector  always 
acts  through  the  aircraft's  center  of  gravity.  This  placement 
means  that  no  pitching  moment  is  caused  by  the  rotation  of 
the  nozzles. 

While  seeing  the  Harrier  in  a  hover  gives  one  the 
impression  that  it  is  defying  gravity,  it  must,  of  course, 
obey  Newton's  laws  of  motion.  This  implies  that  in  all  phases 
of  steady  flight  some  lifting  force  must  balance  the  weight 
of  the  aircraft.  In  the  VTOL/hover  regi.ne,  the  lift  is 
produced  entirely  by  the  thrust  of  the  Pegasus  engine,  the 
nozzles  being  rotated  nearly  perpendicular  to  the  aircraft 
centerline  and  the  thrust  being  directed  vertically  downward. 
For  transition  to  conventional  flight,  the  nozzles,  which 
share  a  common  mechanical  linkage  to  a  control  lever  in 
the  cockpit,  are  slowly  rotated  aft.  The  net  effect  of  this 
rotation  is  to  reduce  the  magnitude  of  the  upward  (lift) 
component  of  the  thrust  vector  while  simultaneously  producing 
a  forward  thrust  component.  The  resulting  forward  velocity 


of  the  aircraft  causes  airflow  over  the  wings  and  therefore 
produces  aerodynamic  lift.  As  the  transition  progresses, 
wing-generated  lift  gradually  replaces  the  upward  component 
of  the  thrust  vector  until,  at  a  speed  of  about  100  knots, 
the  nozzles  are  rotated  fully  aft  and  the  aircraft  is  in 
the  conventional  regime,  The  preceding  sequence  of  events 
is  reversed  for  a  decelerating  transition  from  fully  wing- 
borne  to  fully  jet-borne  fliaht..  It  is  worthwhile  to  note 
here  that  the  Harrier's  most  acute  stability  and  control 
problems  occur  in  the  transition  regime,  from  about  30  to  80 
knots.  The  dynamics  of  this  phenomena  will  be  discussed 
fully  in  the  next  chapter. 

Knowledge  of  the  design  and  performance  of  the  Harrier's 
flight  control  system  is  vitally  important  to  solution  of 
its  stability  and  control  problems.  The  aircraft's  primary 
controls  for  wing-borne  flight  consist  of  conventional  ailerons, 
rudder,  and  stabilator,  wi+h  a  reaction  control  system  (RCS) 
providing  additional  control  authority  in  the  transition 
and  hover  regimes. 

This  reaction  control  system  bleeds  high  pressure  air 
from  the  engine  through  lightweight  ducting  to  special  nozzles 
(called  "puffer"  valves)  located  at  the  nose  tail,  and  each 
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wingtip  of  the  aircraft.  These  shutter  valves  are  mechanically 
linked  to  the  conventional  controls;  a  master  control  valve 
is  geared  to  fully  activate  the  system  when  the  engine  nozzles 
are  rotated  20°  below  the  aircraft's  centerline.  Since  the 
nozzle  angle  is  proportional  to  the  aircraft's  airspeed,  the 
net  effect  of  the  master  control  valve  gearing  is  that  the 
RCS  is  activated  only  in  the  VSTOL-hover  regime  (below  about 
100  KIAS).  The  RCS  extends  the  Harrier's  controllability 
envelope  all  the  way  down  to  zero  airspeed,  long  after  the 
aerodynamic  controls  have  lost  their  effectiveness. 

A  single  channel  limited  authority  stability  augmentation 
system  (SAS)  is  provided  to  dampen  rolling,  pitching,  and 
yawing  motions  while  in  a  hover  or  transition.  The  device, 
which  contains  no  "autopilot1  function,  was  added  to  reduce 
pilot  workload  under  turbulent  and  low  visibility  conditions. 

This  is  accomplished  by  sensing  the  aircraft's  angular  velocities 
about  the  roll,  pitch,  and  yaw  a/.es  and  immediately  initiating 
a  correcting  signal  to  the  control  actuators,  thus  augmenting 
the  pilot's  demands  on  the  aerodynamic  control  surfaces  and/or 
reaction  control  nozzles. 

As  it  is  primarily  a  low  soeed  system,  the  SAS  functions 
almost  exclusively  through  the  reaction  control  system.  It 
should  be  mentioned  that  the  control  authority  of  the  reaction 
control  system  decreases  dramatically  in  high  thrust  demand 
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situations.  The  SAS  is  quite  economical  in  space  and  weight: 
both  the  pitch  and  roll  channels  are  housed  in  a  single  5.5 
pound  unit  that  contains  gyros,  computing,  power  supply,  and 
self-test  devices.  The  yaw  SAS  is  slightly  more  elaborate, 
as  it  utilizes  feedback  of  both  the  aircraft's  yaw  rate  and 
of  its  lateral  acceleration.  In  general,  while  the  addition 
of  stability  augmentation  in  its  present  form  may  definitely 
be  considered  as  a  positive  and  flying-quality-enhancing 
improvement,  there  is  quite  a  lot  of  motivation  toward  opti¬ 
mization  of  the  system. 

Of  additional  interest  from  the  controls  standpoint 
is  a  lateral  acceleration  indication  system  consisting  of 
a  readout  projected  on  the  pilot's  head. -up  display  (HUD) 
and  a  novel  rudder  pedal  shaker  system  that  actually  vigorously 
vibrates  the  rudder  pedal  that  the  pilot  should  depress  to 
eliminate  the  lateral  acceleration.  Sideslip  is  proportional 
to  lateral  acceleration,  and  the  use  of  an  accelerometer  as 
the  sensor  of  a  sideslip-regulation  circuit  for  turn  coordination 
of  conventional  aircraft  is  quite  common.  As  will  be  described 
in  Chapter  3,  the  elimination  of  sideslip  is  of  key  importance 
to  VSTOL  stability  studies.  In  fact,  control  of  this  quantity 
is  so  important  to  the  directional  stability  of  the  Harrier 


1 
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in  transition  that  it  is  fed  buck  not  only  to  the  pilot  through 
the  indication  system,  but  also  to  the  lateral  autostabilizer 
through  an  electronic  filtering/signal  processing  network. 

The  aircraft  that  is  the  subject  of  the  research  effort 
documented  in  this  report  is  an  updated  version  of  the  AV-8A 
now  in  service  in  the  Marine  Corps.  It  is  called 
the  AV-8B  Advanced  Harrier,  and  a  brief  survey  of  its 
development  and  capabilities  is  now  presented. 

Following  the  acquisition  by  McDonnel 1 -Douglas  (1969) 
of  manufacturing  rights  for  the  AV-8A  Harrier,  a  tremendous 
amount  of  research  time  and  money  was  put  into  the  full 
development  of  the  Harrier's  potential.  A  parallel  effort 
was  made  in  the  refinement  of  the  Pegasus  engine  by  Pratt 
and  Whitney  (who  had  similarly  obtained  production  rights 
from  Rolls  Royce).  Many  independent  research  projects 
flourished  on  the  fringes  of  this  central  two-pronged  advance, 
spurred  on  by  publication  of  large  quantities  of  previously 
unavailable  technical  data.  The  AV-8B  Advanced  Harrier  is 
the  improvement  program's  crowning  achievement.  The  linearized 
mathematical  models  developed  by  Calspan  Corporation  of 
Buffalo,  N.Y.,  which  form  the  oasis  for  the  simulation  described 
later  in  this  report,  we»’e  an  offshoot  of  this  engineering 
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renaissance  for  the  Harrier. 


Effective  lift  is  the  AV-8's  most  critical  design  criteria. 
There  are  two  methods  of  increasing  effective  lift:  increase 
lift  or  reduce  weight.  The  AV-8B  makes  use  of  both  approaches. 

To  increase  lift  in  jet-borne  flight,  600  pounds  of  additional 
lift  were  obtained  by  Pratt  and  Whitney  through  a  redesign  of 
the  Pegasus'  engine  inlets.  A  special  "cross  dam"  was  then 
installed  under  the  fuselage  in  order  to  take  better  advantage 
of  interference  effects  during  VTOL.  An  estimated  1200  pounds 
of  additional  lift  is  generated  by  this  device  in  the  VTO 
mode.  Lift  in  the  conventional  flight  regime  has  been  augmented 
by  the  use  of  a  thicker,  "supercritical"  airfoil.  This  new 
wing  also  cuts  down  drag  at  high  speed,  thereby  reducing  fuel 
consumption  -  always  a  concern  in  a  VSTO’  aircraft.  The 
structure  of  the  new  airfoil  is  the  major  weight  reducing 
component  of  the  Advanced  Harrier.  The  entire  wing  is  made 
of  graphite  epoxy  reinforced  at  high  stress  areas  with  titanium 
and  aluminum.  Most  of  the  AV-8B's  other  subsystems  were 
carried  over  from  the  AV-8A  with  minor  modifications.  The 
net  result  of  these  increases  in  effective  lift  is  a  performance 
increase  of  about  100%  (  a  useful  load  of  about  15,000  lbs). 

Like  the  AV-8A,  the  Advanced  Harrier's  operational 
effectiveness  is  best  measured  by  its  versatility.  It  can 


take  over  3000  pounds  of  ordinance  straight  up  and  to  a 
target  fifty  miles  distant;  from  a  short  take-off,  its  per¬ 
formance  characteristics  are  comparable  to  those  of  the  touted 
A-4  Skyhawk.  The  Harrier  has  oroven  itself  capable  of 
operating  from  small  sea-going  platforms,  or  from  dispersed 
austere  sites  ashore.  The  Harrier  is  an  excellent  aircraft; 
it  is,  of  course,  the  only  VSTOL  platform  sufficiently  advanced 
in  development  to  be  used  operationally.  But,  if  any  of  the 
Harrier's  shortcomings  can  be  labeled  "serious",  one  of  those 
must  surely  be  its  transition-regime  stability  characteristics. 
The  next  chapter  examines  this  area  in  depth. 


CHAPTER  3 


STABILITY  CONSIDERATIONS 

A  quick  yet  effective  means  of  forming  a  preliminary 
evaluation  of  an  aircraft's  flight  characteristics  is  to 
note  its  "Prohibited  Maneuvers"  as  tabulated  in  its  flight 
manual.  The  AV-8A  NATOPS  Flight  Manual  [3  ,  R-12]  lists  26 
prohibited  maneuvers,  among  them: 

1)  Out-of-the-wind  vertical  take-off  or  landing  with 
a  wind  speed  greater  than  20  kts. 

2)  In  the  transition  regime: 

a)  Over  15  units  angle  of  attack  above  50  KIAS 

b)  Sideslip  between  30-150  KIAS 

c)  Other  than  gentle  turns  between  90-150  KIAS 

d)  Any  turns  between  30  and  90  KIAS 

3)  Lateral  wind  component  in  excess  of  30  kts. 

4)  Take-off  or  landing  on  paved  surfaces  with  cross- 
wind  greater  than  approximately  5-15  knots  (depending 
on  type  of  landing  and  visibility)  or  on  any  other 
surface  with  crosswind  greater  than  5  kts. 

Thus,  it  is  seen  that  stability  and  control  problems  sub¬ 
stantially  limit  the  Harrier's  usable  flight  envelope, 
especially  in  low  speed  situations. 

The  Harrier's  effectiveness  below  about  80  knots  is  quite 
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limited,  even  when  operating  within  the  boundaries  prescribed 
by  NATOPS.  This  limitation  stems  from  the  unusually  high 
pilot  workload  required  during  the  transition,  especially  in  a 
low  visibility  or  turbulent  environment.  The  present  SAS, 
by  damping  the  unstable,  at  times  erratic  motions  of  the 
aircraft,  very  noticeably  improves  the  Harrier's  handling 
characteristics  in  semi -wing-borne  and  jet-borne  flight. 

But,  as  many  pilots  quickly  attest  to,  an  "...attitude  conmand 
SAS  would  be  desirable  since  this  would  make  a  positive 
orientation  reference  available  to  the  pilot  at  all  times." 

[  4  ,  332]  The  control  system  designed  herein  seeks  to 
answer  this  plea. 

It  is  desirable  at  this  point  to  back-off  from  a  discus¬ 
sion  of  the  aircraft  at  the  limits  of  its  stability/control 
envelope  and  discuss  its  handling  characteristics  in  general 
terms.  First  of  all,  jet  VSTOL  aircraft  operating  in  wing- 
borne  flight  perform  very  much  like  conventional  aircraft 
In  general,  the  Harrier's  flight  characteristics  at  airspeeds 
above  150  knots  are  quite  admirable.  However,  VSTOL  craft 
in  the  jet-borne  regime  may  exhibit  characteristics  that 
appear  quite  peculiar  to  those  familiar  only  with  conventional 
aircraft.  The  Harrier  in  VSTOL  flight,  for  example,  exhibits 
neutral  to  negative  stability  in  pitch  and  roll,  and  a  definite, 
pronounced  instability  in  yaw  due  to  the  intake  momentum  drag 
phenomena.  The  research  effort  described  in  this  report  began 
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as  a  quest  for  a  solution  to  this  particular  control  anomaly 
(along  with  the  related  roll-coupling  and  pitch  sensitivity 
characteristics).  With  this  in  mind,  our  attention  will 
now  focus  on  the  aerodynamic  principles  that  underlie  the 
intake  momentum  drag  effect. 

Air  entering  a  jet  engine  experiences  a  rapid  decrease  in 
velocity  at  the  compressor  section  inlet.  This  reduction 
of  velocity  produces  a  net  change  in  momentum,  which  results 
in  a  "momentum  drag"  force.  At  the  high  engine  power  settings 
required  for  VSTOL  flight,  some  400  Ib/sec  of  air  are  brought 
to  rest  in  the  Pegasus'  intakes.  At  extremely  low  airspeeds 
(below  30  kts),  this  phenomena  has  little  effect,  as  the 
retarding  force  (drag)  caused  by  the  change  in  momentum  is 
quite  small.  At  higher  forward  velocities,  however,  this 
drag  force  is  substantial  (equaling,  according  to  Newton,  the 
product  of  the  mass  of  the  air  and  the  velocity). 

The  air  intakes  of  the  Pegasus  engine  are  located  forward 
of  the  aircraft's  center  of  gravity.  The  geometry  of  the 
situation  implies  that  this  intake  momentum  drag,  in  the  presence 
of  any  angular  deflection  about  the  yaw  axis  (so-called  "sideslip") 
produces  a  destabilizing  moment.  Further,  the  magnitude  of 
this  moment  increases  linearly  with  velocity.  Combating  this 
toroue  is  the  Inherently  stabilizing  yawing  moment  produced  by  the 


aerodynamic  lift  of  the  vertical  tail,  which  increases  as  a 
function  of  the  velocity  squared.  The  net  result  of  these  two 
opposing  moment  mechanisms  is  that  the  Harrier  is  directionally 
very  unstable  between  30  and  90  knots,  but  essentially  stable 
outside  this  velocity  envelope. 

In  the  Harrier,  unless  corrective  control  is  instituted, 
sideslip  disturbances  result  in  coupled  motions  about  all  axes. 

This  is  the  reason  for  the  high  degree  of  emphasis  placed  on 
its  regulation.  The  mechanism  that  underlies  the  aerodynamic 
coupling  associated  with  sideslip  merits  discussion. 

Sideslip  produces  rolling  moments  in  most  aircraft,  especially 
when  swept  wings  are  utilized.  The  rolling  moment  is  generated 
by  a  lift  differential  between  the  two  wings.  This  change  in 
lift  distribution  is  caused,  in  turn,  by  two  mechanisms:  (1) 
the  sideslip  angle  causes  the  airflow  to  pass  nearly  normal  to 
the  leading  edge  of  the  windward  wing,  while  the  leeward  wing  is 
not  only  turned  away  from  the  wind,  but  is  also  shielded  by 
the  fuselage;  and  (2)  the  yaw  rate  produced  by  the  destabilizing 
moment  causes  the  windward  wing's  velocity  relative  to  the  airflow 
to  increase,  and,  conversely,  that  of  the  leeward  wing  to  decrease. 
Thus,  both  mechanisms  augment  the  lift  of  the  wind-ward  wing 
and  reduce  the  leeward  wing's  lift,  thereby  creating  an  unstable 
rolling  moment.  The  righting  characteristic  about  the  roll  axis 
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("dihedral  effect")  that  is  normally  designed  into  an  aircraft 
is  not,  in  most  cases,  of  sufficient  strength  in  the  Harrier 
in  VSTOL  flight  to  combat  the  roll  coupling  instability. 
Furthermore,  since  large  deflections  in  yaw  and  roll  substantially 
reduce  the  control  authority  of  all  flight  controls,  angle  of 
attack  (which  is  analogous  to  pitch  angle,  a  longitudinal  variable) 
is  also  a  coupled  variable  in  this  intake  momentun  drag  situation. 
This  fact  is  particularly  relevant,  considering  the  acute  pitch 
sensitivity  of  the  Harrier  in  transition.  So,  what  began  as  a 
simple  sideslip  disturbance  has,  through  aerodynamic  coupling, 
adversely  affected  handling  qualities  about  all  three  axes. 

The  AV-8A  NATOPS  Flight  Manual  takes  this  analysis  one  step 
further. 

Three  variables  are  of  paramount  importance  to  the  intake 
momentum  drag/roll  coupling/pitch  sensitivity  phenomena: 


Indicated  air  speed,  IAS 

(q) 

Angle  of  attack,  AOA 

(a) 

Sideslip  angle 

(6) 

The  magnitude  of  the  rolling  moment  is  roughly  proportional 
to  the  product  of  the  three  variables  (IAS,  AOA,  $).  [3,  4-7] 
Furthermore,  the  mathematics  are  such  that  if  any  two  of  the 
terms  are  large,  the  rolling  moment  will  be  substantial  even 
for  a  relatively  small  value  of  the  third  term.  Typical  "large" 
values  are: 
q  =  100  KIAS 
a  =  15  units  =  15° 


6  =  30° 
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Note  that  these  numbers  approximately  coincide  with  several 
of  the  limits  mentioned  in  the  "Prohibited  Maneuvers". 


With  respect  to  this  serious  VSTOL  regime  instability, 
the  flight  manual  warns: 


Most  dangerous  of  all,  the  AOA  will  increase  instantly 
with  roll  if  there  is  a  sideslip  angle  present.  This  can 
result  in  an  almost  instantaneous  loss  of  control  with 
very  little  or  no  warning.  A  typical  loss  of  control 
sequence  at  low  IAS  involves  allowing  a  sideslip  to  develop 
which  introduces  a  rolling  moment  which,  if  not  counteracted, 
instantly  increases  AOA  which  increases  the  rolling  moment 
so  that  the  situation  becomes  progressive.  [3,  4-7]. 


The  three  axis  SAS  that  is  presently  installed  in  the  Harrier 
seeks,  through  damping,  to  reduce  pilot  workload  and  lessen  the 
effects  of  disturbances;  the  pilot  would  then  theoretically 
be  able  to  devote  more  attention  to  su^h  concerns  as  sideslip 
control  in  transition.  However,  the  lack  of  an  autopilot  function 
in  the  SAS  puts  a  heavy  burden  on  the  pilot  in  itself,  particularly 
with  respect  to  the  regulation  of  this  yaw/bank/pitch  coupled 
instability  in  turbulent  and/or  low  visibility  conditions. 

Also  the  directional  (sideslip)  controller,  because  of  noise 
problems  associated  with  the  measurement  of  lateral  acceleration, 
can  be  harsh  and  jerking  in  its  operation.  In  addition,  it  is 
extremely  important  to  note  that  simultaneous  application  of 
control  to  more  than  one  reaction  control  system  axis  results 
(In  high  thrust  demand  situations)  in  an  overall  reduction 
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in  the  available  control  authority  about  all  axes.  Recall  that 
the  reaction  control  system  is  the  Harrier's  major  controlling 
agent  in  hover/transition  flight.  This  last  limitation  underlines 
the  fact  that  a  controller  designed  to  combat  this  aircraft's 
characteristics  must  monitor  the  aircraft's  motion  about  all 
axes  and  provide  instantaneous  control  actuation  in  response 
to  disturbances.  The  state  variable  feedback  controller 
described  and  demonstrated  later  in  this  report  does  just  that. 

Having  qualitatively  described  the  Harrier's  configuration 
and  performance,  the  discussion  now  turns  to  a  fairly  rigorous 
mathematical  formulation  of  these  dynamic  characteristics. 
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CHAPTER  4 


THE  EQUATIONS  OF  MOTION 

It  doesn't  take  a  great  deal  of  thought  to  come  to  the 
conclusion  that  the  airplane  in  flight  is  a  very  complex  dynamic 
system.  Its  operatic  requires  a  study  of  six  degrees  of  freedom 
for  gross  motions,  while  the  modes  of  motion  of  its  many  sub-systems 
add  tremendously  to  the  number  of  system  states  required  for  a  model 
of  engineering  precision.  But  a  realistic  model  must,  of  course,  be 
defined,  for  practical  considerations  often  render  it  impossible  for 
the  control  systems  engineer  to  have  any  actual  contact  with  the 
system  he  is  studying.  When  operating  under  this  handicap,  the 
best  he  can  do  is  to  insure  that  he  has  a  thorough  understanding 
of  both  the  mathematical  model  he  is  using  and  of  its  relation 
to  the  physical  world.  With  the  aim  of  providing  a  greater 
degree  of  comprehension  of  the  system  presently  under  study,  this 
chapter  presents  a  derivation  of  the  Harrier  airframe  equations  of 
motion.  The  procedure  used  is  based  largely  on  that  described  in 
Bernard  Etkins',  The  Dynamics  of  Flight,  Editions  I  and  II  [5  and  6], 
an  aircraft  stability  and  control  text  that  has  become  a  standard 
in  the  field. 

In  order  to  provide  a  oerspective  on  the  ensuing  derivation, 
a  suimary  of  the  key  steps  in  the  development  of  the  equations  is 
in  order.  The  aircraft  is  first  treated  as  a  rigid  body  and  the 


-35- 


equations  of  motion  are  derived  assuming  the  reference  axes  are 
fixed  to  its  center  of  gravity.  This  analysis  yields  the  classical 
Euler  equations.  A  discussion  is  then  made  of  the  orientation  of 
the  plane  in  space  and  of  the  choice  of  axes.  The  treatment  of 
spinning  rotors  and  the  modeling  of  the  aircraft  controls  is  then 
considered.  Finally,  small  disturbance  theory  is  used  to  linearize 
and  partially  decouple  the  equations.  It  should  be  emphasized 
that  this  discussion  is  meant  to  be  more  illustrative  than  definitive. 
A  complete  derivation  is  found  in  The  Dynamics  of  Flight. 

A  rigorous  derivation  of  the  general  equations  of  motion  for 

aircraft  begins  with  particle  mechanics,  where  elementary  masses  and 

forces  are  related  via  Newton's  second  Taw  and  summed: 

-  4- 

Z6F  =  £6mV  4.1 

By  applying  the  definition  of  the  center  of  mass  and  carrying  out 

the  summation,  we  arrive  at  the  familiar: 

-*  • 

F  =  m  vc  4.2 

-y 

where  F  is  the  resultant  external  force  applied,  m  is  the  total 
mass,  and  vc  is  the  acceleration  of  the  center  of  mass.  A  similar 
development  for  rotations  equates  the  total  external  moment,  (a,  to 
the  time  derivative  of  angular  momentum,  h  ,  by: 

S  =  ft  4.3 

These  two  vector  equations  (4.2  and  4. 3)  form  the  foundation  for 
the  expanded  equations  of  motion. 
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Rolling  Moment  P  =  Rolling  Velocity  <J>  =  Bank  angle 

Pitching  Moment  Q  =  Pitching  Velocity  0  =  Pitch  angle 

Yawing  Moment  R  =  Yawing  Velocity  <p  =  Azimuth  angl 

X,  Y,  Z  =  Aerodynamic  Force  Components 

U.  V,  W  =  Velocity  Components  (w.r.t.  the  airmass) 

Figure  4.1  Standard  Body-Axis  Nomenclature 


-1  W 

Angle  of  Attack:  =  tan  ^ 

Angle  of  Sideslip:^  =  sin'1  — - — 

vc.g. 

Figure  4.2  Inclination  of  Body  Axes  to  v 


Before  proceeding  with  the  discussion  leading  to  the  Euler 
equations,  a  brief  comment  on  reference  axes  is  appropriate.  First 
of  all,  the  development  presented  here  is  valid  for  any  "body  axis" 
system  --  orthogonal  axes  fixed  in  the  airplane,  with  the  center 
of  gravity  as  the  origin.  Because  of  the  simplifications  they  allow 
later  in  the  derivation,  the  ultimate  choice  of  axes  is  the  "stability 
axis"  system  -  with  the  x-axis  pointing  in  the  direction  of  motion, 
the  y-axis  extending  out  the  right  wing,  and  the  z-axis  pointing 
"downward"  such  that  it  fulfills  the  requirement  of  orthogonality. 
Stability  axes  find  application  primarily  in  the  study  of  linear 
disturbance  theory.  Figures  a.l  and  4.2  define  the  necessary 
nomenclature. 

The  evaluation  of  the  angular  momentum,  h,  is  vitally  important. 

=  £6ii  =  Z(T  x  ^)6m  4-4 

Let  the  aircrafts  angular  velocity  vector  be: 

oj  =  iP  +  jQ  +  kR  4.5 

The  total  velocity  of  a  point  in  a  rotating  rigid  body  is 

-»•  -V  -V  — 

V  =  v^  +  uixr  4.6 

with  r  =  xt  +  yj  +  zli,  from  which 

h  =  Er  +  (v  +  wxr}cm  4.7 

The  application  of  vector  algebra  yields: 

h  =  <fc(x2+y2+z2)5m-i:r(Px+Qy+Rz)6m  4.8 

The  summations  that  occur  in  the  equations  of  the  scalar  components 
of  equation  4.8  are  the  moments  and  products  of  inertia  for  the 
aircraft,  so  we  have,  after  substitution: 
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h' 

"A  -F  -E“ 

V 

X 

h 

= 

-F  B  -D 

Q 

y 

hz 

m  m 

► 

-E  -D  C 

R_ 

where 

[A  B  C]  =  moments  of  inertia  about  x,  y,  and  z  axes,  respectively, 
and  D  =  /yzdm,  E  =  /xzdm,  F  =  /xydm  are  the  so-called  products  of 
inertia. 

The  derivative  of  a  vector  ($  in  a  rotating  reference  frame  is: 
ff  =  j£+wx0  4.10 

where 

V^z 

Applying  this  definition  to  our  vector  equations  of  motion,  in  a 

reference  frame  fixed  to  the  aircraft,  yields: 

6v 

t  =  n>-^  +  mio  x  v  c  4.11 

S  =  |  +  wxii  4.12 

with  the  following  scalar  components: 

Fx  =  m(tl  +  QW  -  R V) 

Fy  =  m(V  +  RU  -  PW)  4  .13 

Fz  =  m(W  +  PV  -  QU) 
and 

L  =  hx  +  Qhz  -  Rhy 

M  =  h  +  Rh  -  Ph,  4.14 

y  x  z 

N  =  hz  +  Phy  -  Qhx 
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where  h  is  defined  by  equation  4.9.  Equations  4.13  and  4.14, 
together  with  4.9,  are  the  Euler  equations  of  motion  for  the 
airplane.  The  physical  interpretation  of  the  variables  of  the 
Euler  equations  will  become  evident  after  the  following  discussion 
concerning  orientation  relative  to  fixed  axes. 

Because  one  of  the  basic  assumptions  made  in  the  foregoing 
analysis  was  that  the  reference  axes  used  were  fixed  to  the  aircraft, 
the  linear  and  angular  position  of  the  aircraft  cannot  be  described 
relative  to  them.  Two  sets  of  differential  equations  produce  the 

l  I  I 

desired  transformation  to  a  fixed  frame  (x  ,  y  ,  z  ): 

•  I 

x  =  UcosOcosiJH-V(sin4)SinecosiJ;-cos<|>sin\J;)+W(cos<t>sinecos4>+sin<{>sin4i) 

y  =  Ucos0sim|^V(sin4>sin0sinijrt-cos<f>cosiJ;)+W(cos4>sinesin4>-sin<j>cosH;)  4.15 

z'  =  -Usin0+Vsin<}>cos0+Wcos<{>cos0 

and,  to  describe  the  angular  orientation: 

•  • 

P  =  <f>  -  if  sin  0 

Q  =  0  cos  <J>  +  ip  cos  0  sin  <p  4.16 

R  =  ip  cos  0  cos  <p  -  0  sin  4> 

therefore, 

0  =  Q  cos  <f>  -  R  sin  <j> 

<j>  =  P  +  Q  sin  <b  tan  0  +  R  cos  4>  tan  0  4.17 

ip  =  (Q  sin  4>  +  R  cos  $)  sec  0 

It  is  obvious  that  the  integration  of  equations  4.15  and  4.17  represents 
a  formidable  problem  except  in  simplified  special  cases.  Fortunately, 
linearization  by  small  disturbance  theory  simplifies  them.  It  should 
be  mentioned  that  the  orientation  of  an  airplane  relative  to  another 


axis  system  may  be  uniquely  defined  by  an  ordered  series  of  three 
rotations  -  first  in  yaw  to  angle  ifj,  then  in  pitch  to  angle  6, 
and  finally  in  roll  to  angle  <J>.  The  specific  ordered  set  of 
rotations  just  described  defines  the  so  called  Euler  angles. 

The  next  step  in  the  development  of  the  equations  is 
to  account  for  the  "rotary  derivatives"  -  the  angular  momentum 
added  to  the  total  by  any  spinning  subsystems.  This  quantity 
would  normally  appear  as  an  additive  term  (li1)  to  the  r.h.s.  of 
equation  4.8,  but,  because  of  the  Harrier's  contrarotating  turbine 
sections,  it  is  not  important  to  this  discussion. 


The  next  subsystem  group  to  be  modeled  is  that  of  the  flight 
controls.  A  detailed  discussion  of  this  area  is  beyond  the  scope 
of  this  report.  In  general,  each  of  the  three  control  systems  (rudder, 
aileron  and  elevator)  is  assumed  to  be  a  rigid,  one-degree-of-freedom 
linkage.  LaGrange's  equation  of  motion  in  a  moving  reference  frame 
is  then  applied: 


d_  dT  dT  = 
dt  3qk  "  3qk 


4.18 


where 

T  =  kinetic  energy  of  the  system 
W  =  work  done  on  the  system  by  external  forces 
F  =  generalized  force 
qk  =  generalized  coordinate 
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Application  of  this  equation  yields  the  three  control  equations: 

He  +  Fe  =  Ve  *  "eVc,  *  Pex  <PR-'I>  (alevator> 

Hr  ♦  Fr  =  Ir5r  -  VraCy  '  Prx  (PQtP)  *  Prz  (R0'P)  (rudder>  4-19 

2H  +  F  =  r  *6*  +  2P  (RQ  +  P)  (aileron) 
a  a  a  a  ay 

where 

a  =  accelerated  of  center  of  mass 
c 

e  =  mass  eccentricity  of  control  surface 
F  =  generalized  control  force 
H  =  control  hinge  moment 

I  =  control  effective  moment  of  inertia 

M  =  mass  of  control  surface 

P  =  control  product  of  inertia 

6  =  control  angle 

The  development  of  the  control  equations  represents  the  last  step 
in  the  modeling  of  the  aircraft  from  a  stability  and  control  point  of  view. 
However,  the  complexity  of  the  equations  renders  them  nearly  unusable 
in  their  present  form. 

One  minor  addition  needs  to  be  added  to  the  present  equations. 

The  total  external  force  applied  to  the  aircraft  is  the  sum  of 

aerodynamic  (including  propulsive)  and  gravitational  forces.  The 

quantities  appearing  on  the  l.h.s.  of  equations  4.13  are,  therefore, 

F  =  X  -  mg  sin  0 
x  3 

Fy  =  Y  +  mg  cos  6  sin  <p  4.20 

Fz  =  Z  +  mg  cos  0  cos  <p 


-42- 


The  moment  equations,  being  purely  aerodynamic,  need  not  be  altered. 

Small  disturbance  theory  will  now  be  used  to  linearize  the 
equations.  The  validity  of  the  approach  rests  on  two  facts: 

(1)  it  is  often  the  case  that  the  aerodynamic  response  is  a  nearly 
linear  function  of  the  disturbance,  and  (2)  flight  involving  large 
external  disturbances  may  be  sustained  with  relatively  small  variations 
in  linear  and  angular  velocity  of  the  aircraft.  Small  disturbance 
theory  necessitates  the  following  assumptions: 

(a)  all  disturbances  and  derivatives  all  small,  and  their 
products  and  squares  can  therefore,  be  considered 
negligible; 

(b)  stability  axes  are  used; 

(c)  reference  flight  condition  symmetric  and  without 
angular  velocity; 

(d)  small  angle  approximations  may  be  used; 

(e)  rotary  derivatives  are  negligible;  and 

(f)  all  controls  are  dynamically  balanced. 

The  application  of  small  disturbance  theory  to  linearize  the 
equations  involves  several  steps.  Initial  application  to  equations 
4.13,  4.14,  4.16,  4.17,  and  4.19  yields: 

XQ  +  AX  -  mg(sin0o  +  0cos0o)  =  mu 

Yq  +  AY  +  mg4icos0o  =  m('vf  +  uQr)  4  .21 

Z^  +  AZ  +  mg(cos0„  -  0sin0  1  =  m(u>  -  u  a) 

O  O  0  0 


I 


•  • 

LQ  +  AL  =  Ap  -  Er 

M_+  AM  =  Bq  4.22 

0 

NQ  +  AN  =  -Ep  +  Cr 


+  AH  + 

F  +  AF  = 

In 

e 

eo  e 

e 

+  AH  + 

F.  +  AF  = 

I  a 

4  .23 

r 

ro 

r 

+  2AH 

+  F  +  AF 

=  i  i' 

>  a 

aO  3 

a 

0  =  q 

4>  =  p  +  r  tan  0Q  P  =  <J>  -  4>  sin  dQ  4.24 

\p  =  r  sec  0Q 

If  the  disturbance  variables  in  equations  4.21  -  4.24  are  set  equal 
to  zero,  the  equations  then  apply  to  the  reference  condition.  By 
subtracting  these  reference  values  from  the  equations,  we  can 
eliminate  all  initial  forces  and  moments: 


AX-mg0cos0Q=  mu 

AY+mg<j>cos0o=  mv  +  muQr  4 .25 

AZ-mg0sin0Q  =  mu  -  muQq 

AL  =  Ap  -  Er 

AM  =  Bq  4 . 26 

AN  =  -Ep  +  Cr 
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I 

I 

AHe  +  AFe  =  !pN 

T 

AHr  +  AFr  =  Ir*C  4.27 

2AH  +  AF  =  l’i 

a  a  a 

e  =  q 

<p  =  p  +  r  tan  6  p  =  <f>  -  sin  0  4  .24 

(unchanged) 

ip  =  r  sec  9q 

The  derivation  up  to  this  point  has  been  general  in  scope, 
i.e.,  is  applicable  to  both  conventional  and  VSTOL  aircraft. 

It  is  in  the  treatment  of  the  aerodynamic  force  and  moment 
perturbations  (the  terms  with  a  "A"  prefix)  that  the  differences 
show.  This  is  also  the  point  of  major  problems  from  the  VSTOL 
modeling  point  of  view,  as  the  literature  essentially  omits 
this  process  from  their  discussion.  The  actual  method  used  was 
introduced  by  Bryan  in  1911.  It  is  based  on  the  assumption 
that  the  aerodynamic  perturbations  can  be  constructed  from 
the  instantaneous  values  of  the  disturbance  velocities,  the 
control  deflections,  and  their  derivatives  via  a  Taylor 
series  expansion.  The  resulting  series  is  linearized  by 
neglecting  all  higher  order  terms.  One  such  "linear  air 
reaction"  takes  the  form 

M  =  Auu  +  Auu  +  ...A66  +  A-6  4.28 

where  A,,  =  (-?-)  reference  condition 

U  oil  0 
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The  derivatives  "A"  are  called  the  stability  derivatives 
of  the  aircraft.  Engineering  experience  and  judgement  are 
used  to  determine  which  stability  derivatives  are  important 
enough  to  be  included  in  the  air  reaction.  This  decision, 
in  turn,  affects  the  final  form  of  the  equations.  Further 
discussion  of  this  decision  process  is  beyond  the  scope 
of  this  text. 


After  careful  thought  and  consideration  of  the  literature, 
it  was  concluded  that  the  VSTOL  linear  air  reactions  take  the 
following  form  for  the  Harrier: 


M  -  Xuu  ♦  Xww  ♦  Xqq  ♦  X^  0j  ♦  X{t«t  ♦  X^ 

"  '  v *  V  *  V  +  \SA  *  Y{rSr 


4Z  =  Zull  ♦  Zww  +  Zqq  +  Z^Oj  *  Z^ST  +  Z^Sj 


4L  *  V  *  Lpp  +  Lrr  +  L6.SA  *  L6.6R 

A  K 

AM  =  M  u  +  M  w  +  Mo  +  Mfl  6 .  +  M.  ST  +  M,  6C  4 .30 

u  w  q^  6j  J  Oj  I  E 

AN  =  Nyv  +  N  p  +  Npr  +  Ng  +  Ng  &r 
v  A  r 

and  I.  ■  I-  *  I.  *  0 
are 
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Substitution  of  the  air  reactions  into  equations  4.24-4.27 
yield  equations  similar  in  form  (but  not  content)  to  those 
in  the  literature: 


CmK  -  xu]  u  -  V  *  <m9  C0S60  -  9  -  *6.6,  -  *s  «T  -  *s  «e  -  o 

J  T  e 

"Zuut(m^  ’  zw)w  ~  K"V  Zq^  ‘  m9  sine0]  9  -  ze.ej  -  z6  6t  ~h  6e 

j  T  e 


-Huu-Hwwt(B^7  -  M,  7E>e  -M6.9j  -  MSt«t  -  Wje6e  ■  o 


q  -*e  0 


^  “  V  v  "  YpP  +  (mu0-Yr)r-(mg  cosfjtf  -Y^  <5a-Y6  6r  =  0 

a  r 

-Lvv+<a£  -  Lp)p-(E^  ♦  Lr)r-Ls>Sa-L{rSr  =  0 

"Nvv  -  tilt  *  "p>P  *  -  Nr)r  -  N{  «  -  N  &  =  0 

a  r 

<t>  =  p  +  tan  9Qr 
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OUTPUT  EQUATIONS 


ot  =  tan  1  J  =  tan  [7^] 

0 

-1  V  -1 

3  =  tan  y\  s  tan  r_ ^ _ 1 

w0+u 

ij/  =  r  sec  6Q 


4.33 


Since  9=5°  for  the  model,  the  second  term  of  the  r.h.s. 

0 

of  the  last  equation  of  equation  set  4.32  may  be  neglected. 
Algebraic  manipulation  is  all  that  is  needed  to  transform 
equations  4.31  and  4.32  into  their  final,  state  variable 
matrix  form.  The  results  are  shown  in  Figures  4.3  and  4.4. 

But  equations  4.31  and  4.32  are  of  little  use  without  experi¬ 
mental  knowledge  of  the  values  of  their  coefficients.  Hence, 
the  next  chapter  discusses  the  source  from  which  the  AV-8B's 
stability  derivatives  were  obtained,  and  describes  the  methodology 
used  to  determine  them  experimentally. 


CHAPTER  5 

ORIGIN  OF  THE  MODEL 

The  simulation  described  later  in  this  report  is  based  upon  a 
linear,  non-constant  coefficient  mathematical  model  developed  by 
Cal  span  Corporation  of  Buffalo,  New  York  in  the  summer  of  1977. 

[  7]  Cal  span  actually  laid  the  groundwork  for  the  model  much 
earlier,  publishing  a  feasibility  study  which  contained  a  generic 
VSTOL  model  in  July  of  1976.  [  8]  The  genesis  of  the  linear  AV-8B 
model  as  published  in  the  final  report,  Calspan  TM  No.  98,  is  in¬ 
teresting  in  two  major  respects. 

The  first  is  that  the  original  motivation  for  developing  a 
linear  math  model  was  to  determine  the  feedback  gains  for  the 
variable  stability  mechanism  of  Calspan's  X-22A  in-flight  simula¬ 
tion  aircraft.  In  fact,  hardware  limitations  of  the  X-22A  avionics 
were  the  primary  reason  for  linearization  of  the  AV-8B  equations 
of  motion.  The  other  particularly  notable  fact  concerning  the 
origin  of  the  model  is  that  it  is,  essentially,  a  simulation  of  a 
simulation.  Simulation  of  the  actual  Advanced  Harrier  was  not 
possible  because  at  the  time  of  the  TM  No.  98  report,  the  AV-8B 
was  not  in  existence.  A  study  of  the  methodology  used  in  the 
TM  No.  98  model  development  provides  needed  insight  as  to  the 
fidelity  and  application  of  the  linearized  equations  of  motion 
for  the  Harrier  used  in  this  work. 

The  only  model  of  the  AV-8B  Harrier  in  existence  at  the  time 
of  the  Calspan  report  was  a  nonlinear  table-look-up  digital 
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computer  simulation  at  McDonnell -Douglas'  (MCAIR)  St.  Louis 
facility.  Consequently,  Calspan,  in  cooperation  with  NADC, 

Warminster,  requested  that  two  representative  approach  (transition-* 
hover)  trajectories  be  "flown"  on  the  MCAIR  simulation.  The  air¬ 
craft  was  to  be  kept  at  trim  throughout  the  approach.  MCAIR  com¬ 
plied,  defining  a  3°  (flight  path  angle  or  "glide  slope")  ap¬ 
proach  from  105  Kt,  and  a  5°  approach  from  65  Kt.  It  was  decided 
that  the  angle  of  attack  would  be  held  constant  at  8°  and  that  a 
"one-step"  deceleration  would  be  used,  that  is,  nozzles 
rotated  from  full  back  to  full  down  instantaneously.  In  all,  these 
two  trajectories  may  be  considered  representative  of  actual  flight 
procedures  based  on  the  NATOPS  Flight  Manual.  [  3] 

The  mathematical  model  was  to  be  developed  in  terms  of  these  task- 
related  reference  trajectories  ,  with  linearized,  but  non-con¬ 
stant  stability  and  control  derivatives  being  calculated  about 
the  "trim"  conditions.  Small  perturbation  theory  and  the  use 
of  stability  axes  as  the  reference  frame  are  necessary  to  employ  tnis 
procedure.  The  actual  calculation  of  the  derivatives  was  conducted 
as  follows. 

A  computer  subroutine  was  written  to  interact  with  the  MCAIR 
simulation.  At  various  points  along  each  reference  trajectory, 
the  subprogram  first  perturbed  one  of  the  state  variables  (u,  v, 
w,  etc.),  and  then  measured  the  corresponding  changes  in  the  other 
aerodynamic  and  thrust  forces  and  moments  (the  other  state 
variables).  These  force  and  moment  deflections  were  then  plotted 


versus  the  perturbation,  and  the  appropriate  linearized  stability 
derivative  determined  by  statistical  data  fit  in  the  least  squares 
sense.  Since  this  entire  analysis  is  based  on  perturbation 
studies  at  a  single  velocity  along  the  reference  trajectory,  the 
process  yields  a  discrete  series  of  stability  derivatives,  one  for 
each  velocity  tested.  Calspan  choose  0,  30,  50,  65,  80,  and  105 
knots,  and  TM  No.  98  contains  tabulated  values  of  every  stability 
and  control  derivative  at  each  of  these  six  velocities- [  7  , 
pp  9-10]  A  summary  of  these  results  is  contained  in  Appendix  A  . 

In  terms  of  the  equations  of  motion,  the  stability  derivatives 
represent  the  coefficients  of  the  characteristic  equation.  As 
shown  in  the  last  chapter,  a  linear  model  can  be  represented  as  a 
set  of  simultaneous  first  order  differential  equations  of  the  form 

^=£X+£G  51 

where 

A  =  n  x  n  Coefficient  Matrix 

JB,  =  n  x  r  Input  Matrix 

l  =  n  x  1  State  Vector 

u  =  r  x  1  Input  Vector 

The  stability  and  control  derivatives  are  then  the  elements  of  the  a 
and  B  matrices  respectively. 

Thus,  the  airframe  model  of  the  AV-8B  Harrier  has  been  com¬ 
pletely  developed  and  defined.  It  is  seen  to  consist  of  two  fourth 


order  systems  of  the  form  of  equation  5.1.  One  equation  describes 
the  lateral /directional  modes  of  motion,  while  the  other  treats 


the  longitudinal  characteristics.  The  coefficients  of  the  A 
and  B  matrices  have  been  calculated  by  Calspan  for  six  distinct 
points  in  the  transition/hover  regime,  and  have  been  shown 
to  be  functions  of  aircraft  velocity.  Having  derived  and 
defined  a  Harrier  airframe  model,  the  discussion  can  now  focus 
on  the  theory  of  stability  augmentation  controllers. 
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CHAPTER  6 


PRESENT  STABILITY  AUGMENTATION 

It  is  appropriate  to  preface  an  analysis  of  the  way  the 
present  Harrier's  Stability  Augmentation  System  (SAS)  attempts 
to  solve  the  instability  problem  in  transition  flight  with  a 
few  comments  on  the  development  of  the  SAS  block  diagram  used 
in  this  research.  In  the  first  place,  no  definitive  model  of 
the  AV-8B  SAS  was  found  in  the  literature.  The  most  thorough 
treatment  of  the  AV-8A  SAS  is  contained  in  Reference  8  as  an 
appendix.  [8,  pp  92-5]  NADC,  Warminster  responded  to  the 
writer's  request  for  more  information  by  providing  a  sketch  of 
their  best  estimate  of  the  SAS  structure.  The  block  diagrams 
of  the  AV-8B  SAS  shown  in  Figures  6.1  and  5.2  were  constructed 
using  information  from  various  sources,  judgements  being  made 
based  on  engineering  experience  in  cases  of  conflicting  data. 

This  updating  of  the  model  was  a  dynamic  and  on-going  process 
throughout  the  research.  The  final  models  arrived  at  and  shown 
here  are  believed  to  accurately  represent  the  AV-8B  SAS. 

Inspection  of  Figures  6.1  and  6.2  reveals  that  all  channels 
of  the  SAS  exhibit  a  saturation  nonlinearity.  These  nonlinearities 
arise  from  the  control  authority  limits  imposed  on  the  SAS. 


Nonlinear  system  behavior  cannot  be  compensated  by  classical 
controls  techniques.  In  addition,  it  can  be  seen  that  the  SAS 
model,  by  requiring  complex  signal  filtering  and  processing 
in  the  feedback  loops,  significantly  raises  the  order  of  the 
system.  Also,  except  for  the  lateral  acceleration  feedback 
loop  (which  will  be  discussed  below),  angular  rates  are  the 
only  system  variables  fed  back.  Feedback  of  the  first  time 
derivative  of  the  controlled  quantity  usually  effects  relative 
stability,  serving  to  increase  the  damping  of  the  motion.  Thus, 
Harrier  pilots  commonly  refer  to  the  SAS  as  a  "control  damper". 

In  the  literature,  angular  rate  feedback  in  aircraft  controllers 
without  accompanying  feedback  of  angular  position  is  virtually 
nonexistent. 

The  feedback  of  single  states  to  single  -'nputs  with  accompanying 
series  compensation  networks  in  a  multi-input,  multi-output 
system  is  an  attempt  to  apply  something  familiar,  classical 
transfer  function  compensation,  to  something  unfamiliar,  the 
multi-input,  multi-output  plants.  The  present  Harrier  SAS 
represents  just  such  an  attempt.  It  seems  safe  to  conclude 
that  this  approach  can  provide  neither  a  systematic  design  of, 
nor  an  optimal  configuration  of  total  system  compensation. 

Because  the  Harrier  lateral  and  longitudinal  airframe  models 
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consist  of  coupled  differential  equations,  and  because  both 
systems  of  equations  are  of  the  multi-input,  multi-output  type, 
conventional  transfer  function  control  design  is  largely  in¬ 
effective.  This  is  not  to  say  that  it  is  impossible  to  determine 
the  input/output  relation  (transfer  function)  between  any 
given  input  and  output.  With  a  knowledge  of  the  composition 
of  the  £  and  B  matrices,  it  is  relatively  easy,  with  the  use  of 
modern  control  theory,  to  perform  the  necessary  matrix  mani¬ 
pulations  for  transfer  function  determination.  The  problem 
that  arises  is  that  compensation  of  the  determined  input-output 
network  usually  has  a  detrimental  effect  on  the  performance  of 
the  other  input-output  relationships.  In  the  interest  of 
completeness,  this  method  of  SAS  design  was  attempted  with  the 
Harrier  airframe  equations  using  the  powerful  Thayer  method 
of  classical  transfer  function  compensation.  As  might  be 
expected,  the  attempt  was  no  more  successful  than  the  present 
SAS. 


From  a  performance  standpoint,  the  present  SAS  design  does 
improve  the  Harrier's  handling  qualities  to  a  limited  extent. 
From  a  state-of-the-art  control  systems  point  of  view,  however, 
the  system  is  quite  rudimentary,  certainly  yields  suboptimal 
results  (the  aircraft  is  still  unstable),  and  is  therefore 
inadequa te. 


CHAPTER  7 


STATE  FEEDBACK  CONTROL 

It  seems  intuitively  correct  that  the  more  knowledge  the 
controller  has  as  to  the  state  of  the  system  being  controlled, 
the  more  accurate  the  controlling  signal  produced.  Similarly, 
if  the  controller  has  the  ability  to  affect  a  greater  number 
of  control  inputs,  it  will  be  better  able  to  bring  the  controlled 
plant  to  the  desired  output  state.  By  measuring  every  system 
output  state,  and  then  calculating  an  actuating  signal  for 
each  control  input  of  the  system  based  on  a  specified  weighted 
sum  of  the  output  states,  the  state  feedback  controller  generates 
a  superior  solution  to  the  situation  discussed  in  the  previous 
chapter.  This  chapter  documents  the  development  of  the  digitally- 
based,  state  feedback  controller  designed  a:.  a  stability  augmentation 
device  for  the  AV-8B  Harrier.  A  representative  drawing  of  the 
system  is  presented  in  Figure  7.1. 

The  use  of  state  feedback  control  eliminates,  or  greatly 
reduces  all  the  problems  inherent  to  the  design  of  the  present 
SAS.  To  review,  the  deficiencies  of  the  present  SAS  are: 

(1)  failure  to  account  for  nonlinear  effects; 

(2)  requirement  of  complex  analog  signal  processing  in 
the  feedback  loop; 

(3)  greatly  increased  system  order; 

(4)  has  no  significant  effect  on  root  locations;  and 

(5)  does  not  account  for  all  the  states,  nor  for  coupling 
effects. 
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Figure  7.1  Representative  Model:  Digital  State  Feedback 


Figure  7.2  State  Feedback  Block  Diagram 


These  problems  are  not  encountered  in  the  state  feedback 
controller.  In  addition,  the  following  benefits  are  realized 
using  state  feedback: 

(1)  Because  the  controller  reacts  so  efficiently  to  counter¬ 
act  disturbances  of  normal  magnitude,  saturation  is 
seldom  encountered,  and  the  system  is  thus  effectively 

1  inear. 

(2)  Output  variables  are  multiplied  by  a  pure  gain  and  summed, 
thus  no  extensive  signal  processing  is  required. 

(3)  System  order  is  unaffected  by  state  feedback  control. 

(4)  Provided  no  limitations  are  put  on  feedback  gains, 
it  is  theoretically  possible  to  relocate  closed  loop 
poles  anywhere  in  the  s-plane. 

(5)  The  calculation  of  the  state  feedback  gain  matrix,  kT, 
accounts  for  coupling  effects  between  all  system  states 
and  inputs. 

Basic  modern  control  theory  states  that  the  characteristic 
differential  equation  for  a  system  can  be  determined  by  setting 
the  determinant  of  the  [sNA]  matrix  equal  to  zero.  The  roots 
of  the  resulting  algebraic  equation  in  the  complex  variable  "s" 
are  then  equal  to  the  roots  of  the  system  itself,  and  therefore 
define  the  system's  natural  modes  of  motion.  When  this  analysis 
is  applied  to  the  Harrier  airframe  equations  (see  Appendix  B), 
it  is  found  that  at  least  one  root  of  the  lateral  equation  and 
one  root  of  the  longitudinal  equation  has  a  positive  real  part. 
This  is  to  be  expected  as  the  flight  manual  contains  numerous 
warnings  of  transition  instabilities.  The  function  of  a  stability 
augmentation  system  (SAS)  then,  in  terms  of  the  roots  of  the 
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characteristic  equation,  is  to  "pull"  the  loci  of  root  positions 
(or  "root  loci")  further  to  the  left  and,  if  possible,  into 
the  stable  left  half  plane.  The  design  criteria  used  to 
determine  the  SAS  state  feedback  gains  for  this  study  was 
the  placement  of  the  system  closed  loop  poles.  The  following 
discussion  is  a  brief  summary  of  the  technique  as  presented 
in  Reference  9.  [9,  pp.  308-314] 

It  can  be  shown  that  for  a  controllable  open  loop  system 
of  the  form 


x  =  Ax  +  Bu 

any  desired  closed  loop  poles  may  be  achieved  by  a  constant  st..r 
feedback  matrix  kT.  It  can  also  be  shown  that  the  numerator 
dynamics  of  the  closed  loop  system  are  unaffected  by  state 
feedback.  Furthermore,  the  closed  loop  eigenvalues  must  satisfy 
the  equation  (refer  to  Figure  7.2): 

A 1 ( A)  =  | X i I  -  A  +  BkT|  =  0  7.2 

for  each  desired  root,  A^ .  It  is  possible  to  transform  this 
equation  into  an  equivalent  form  (see  "Nomenclature"  section  for 
explanation  of  symbology): 
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4'(Xf)  .  0  -  l(X.I„-A)[In  +  <A(In-A)-'BkT|| 

"  +  (Vn-*>''i*kT| 

*  *  <t(x)BkT| 

since  A^.)  =  lXiIn"Al 

=  A(A.)|Ir  +  ?4>(Ai )  BT|  7  .4 

Here  4>(A^)  =  (A-In-A)  and  A(A^)  is  the  open  loop  characteristic 
equation.  The  kT  that  drives  the  determinant  specified  by  either 
equation  7.2  or  7.3  to  zero  for  each  X.  is  the  solution  to  the  pole 
placement  problem. 

Elementary  matrix  theory  insures  a  zero  determinant  if  any 
row  or  column  is  zero.  Associate  Professor  E.  E.  Mitchell,  the 
principal  faculty  advisor  for  this  project,  programmed  an 
algorithm  to  force  the  columns  of  equation  7.4  to  zero.  A 
listing  of  this  program  is  included  as  an  appendix. 

After  suitable  manipulation  of  the  equations,  it  is  possible 
to  determine  that  there  are  nr  unique  solutions  to  the  feedback 
gains  problem.  This  is  fortunate  as  the  algorithm  proved 
susceptible  to  computer  round-off  error  noise  in  this 
application.  (In  general,  the  round-off  error  is  more 
prominant  when  the  plant  is  a  large,  unstable  system).  To 
combat  the  round-off  problem,  it  is  good  practice  to  use 
double  precision  variables  and  to  test  for  a  zero  determinant 
for  each  solution. 
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The  final  design  procedure  was  as  follows.  First,  the  gain¬ 
calculating  program  (see  Appendix  D)  was  employed  to  study  the 
magnitudes  of  the  gains  required  for  various  pole  placement 
areas.  While  infinite  gain  is  theoretically  possible  in  a 
digital  controller,  extremely  high  gains  are  still  undesirable. 
Time  responses  for  the  various  gain  matrices  were  studied  on 
the  digital  computer  as  a  part  of  this  effort.  As  a  result, 
the  following  system  pole  locations  were  somewhat  subjectively 
chosen: 

Longitudinal:  -2,  -2.2,  -2.4,  -3 

Lateral:  -3,  -3.2,  -3.5,  -4 

Not  only  did  these  roots  give  satisfactory  system  time  responses, 
but,  for  Vcq  =  30  kt,  they  were  small  enough  to  be  magnitude 
scaled  to  run  on  the  analog  section  of  a  hybrid  computer.  There 
are  two  major  reasons  for  the  desirability  of  selecting  roots 
on  the  real  axis:  (1)  It  was  thought  that  roundoff  errors 
would  make  the  pole  placement  a  "rough"  process,  and  so  a  high 
damping  ratio  was  desired;  and  (2)  the  gain-calculating  algorithm 

proved  especially  unstable  for  complex  pole  placement.  Like 

the  stability  derivatives,  the  state  feedback  gains  are  also 
velocity  dependent.  A  table-look-up  routine  was  thus  utilized 
to  implement  variable  gains  as  functions  of  velocity. 

On  paper,  state  feedback  control  provides  a  plethora  of 
advantages  including:,  (1)  positive  placement  of  poles  in  the 
s-plane;  (2)  "total  system"  compensation  accounting  for  coupling 
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effects;  and  (3)  no  real  increase  in  model  complexity.  The 
sunmary  of  simulation  results  leaves  little  doubt  as  to  the 
validity  of  the  technique. 
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CHAPTER  8 

COMPUTER  SIMULATION 

The  final  model  of  the  Harrier  airframe  consists  of  two 
fourth  order,  linear  differential  equations  of  the  form: 

x  =  Ax  +  Bu  8.1 

with  the  elements  of  the  A  and  B  matrices  being  the  velocity- 
dependent  stability  and  control  derivatives,  respectively. 

The  first  step  toward  a  computer  simulation  of  this  set  of 
equations  was  the  development  of  a  subroutine  to  calculate 
and  read  the  stability  derivatives  at  any  velocity  (see  Append: x 
E).  Since  the  values  of  the  derivatives  were  known  only  for 
six  distinct  velocities,  some  sort  of  interpolation  scheme 
was  required.  The  method  chosen  was  a  polynomial  fit  in  the 
least  squares  sense.  The  coefficient-calculating  subprogram, 
MATCOF,  reads  the  coefficients  of  the  least-squares  polynomials, 
and  then  uses  them  to  calculate  the  stability  and  control 
derivatives.  The  derivatives  are  returned  in  the  arrays  ALON, 
BLON,  ALAT,  and  BLAT.  This  subprogram  was  used,  without 
modification,  on  both  the  large  scale  digital  computer  simulation 
and  the  hybrid  simulation. 

The  initial  testbed  for  the  mathematical  model  was  the 
Naval  Academy's  Time-Sharing  System  (NATS).  A  FORTRAN  digital 
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simulation  utility  program,  DIGISIM,  was  used  as  the  basis 
for  the  simulation.  DIGISIM,  which  was  also  programmed  by 
Associate  Professor  Mitchell,  uses  a  predictor-corrector  algorithm 
to  handle  the  integrations.  The  many  user-defined  options 
available  in  DIGISIM  make  it  a  particularly  powerful  simulation 
tool.  These  features  are  described  in  Appendix  G  of  this  report. 
Using  an  existing  utility  program  has  the  additional  advantage 
of  some  assurance  of  the  stability  of  the  integration  routine 
itself. 

The  DIGISIM-based  digital  simulation,  named  AHEAD  -  Advanced 
Harrier  Electronic  Augmentation  Device  (see  Appendix  H),  was 
a  tremendous  aid  to  the  development  of  both  the  model  of  the 
present  SAS  and  the  analog/hybrid  simulation.  It  proved 
particularly  helpful  in  the  determination  of  scaling  constants 
for  the  analog  implementation  of  the  equations  of  motion.  It 
was  also  used  to  check  the  accuracy  of  the  analog  model.  There 
is  every  reason  to  believe  that  AHEAD  is  an  excellent  computer 
model  of  the  AV-8B,  its  present  SAS,  and  the  newly-designed 
digital  state  feedback  controller.  The  major  drawback  to  its 
use  as  a  controller-design  tool  is  that  the  simulation  consumes 
approximately  ninety  seconds  of  NATS  C.P.U.  time  for  ten  seconds 
of  run  time. 


Because  it  possesses  (in  its  analog  section)  the  ability 
to  integrate  in  real  time  coupled  with  no  loss  of  the  con¬ 
venience  of  digital  program  control  and  input/output,  the  hybrid 
implementation  of  the  Harrier's  equations  of  motion  represents 
the  ultimate  tool  in  controller  design.  Because  of  the 
complexity  of  the  system,  this  implementation  was  a  long  and 
painstaking  process.  Noise  in  the  analog  section  is  by  far 
the  most  formidable  problem  to  be  overcome  when  modeling  a 
large,  unstable  system.  Persistance  paid  off,  however,  and  the 
final  version  of  the  hybrid  simulation  yielded  responses  that 
compared  closely  to  those  of  AHEAD  for  simulation  times  up  to 
four  seconds  for  the  lateral  and  longitudinal  equations 
(maximum  error  recorded  was  approximately  29%). 

Noise  is  the  limiting  factor  as  far  as  accuracy  is  concerned, 
especially  in  the  longitudinal  simulation.  Two  items  should 
be  noted  regarding  simulation  accuracy:  (1)  the  noise  had 
little  effect  on  the  position  variables,  probably  due  to  the 
filtering  effect  of  the  pure  integration;  and  (2)  a  maximum 
error  of  29%  on  one  state  of  an  unstable  fourth  order  system 
is  actually  quite  small.  Thus,  the  hybrid  simulation  is  also 
judged  to  yield  an  excellent  representation  of  the  Harrier's 
flight  in  the  VSTOL  regime. 


The  one  really  major  problem  encountered  in  the  analog/ 
hybrid  simulation  involved  the  lateral  acceleration  feedback 
loop  of  the  present  SAS  model.  In  order  to  generate  the 
lateral  acceleration  on  the  analog  model,  the  derivatives  of 
three  of  the  states  were  required.  The  derivative  is  available 
by  summing  the  inputs  to  the  particu’ar  state's  integrator. 

This  was  tried,  but  the  resulting  lateral  acceleration  signal 
was  so  noisy  that  its  feedback  actually  had  a  destabilizing 
effect  on  system  response.  The  circuit  was  thereafter  removed 
from  the  hybrid  model.  It  is  interesting  to  note  that  the 
actual  lateral  accelerometer  used  in  the  Harrier  is  very  sensitive, 
having  a  threshold  acceleration  of  0.06  6.  As  such,  this  signal 
is  highly  sensitive  to  spurious  system  motions.  In  addition, 
the  literature  suggests  that  the  feedback  of  lateral  acceleration 
is  usually  accompanied  by  noise  problem.  [1,  p.  145  and  2,  p.  17], 
This  report  supports  that  belief. 

The  final  hybrid  model  represents  a  truly  versatile  tool 
for  the  evaluation  of  control  system  performance.  The  initial 
analog  parameters  (including  potentiometer  settings)  are  set 
by  a  digital  program  that  shares  its  name  with  its  NATS  brother, 
AHEAD  (see  Appendix  J).  This  program  automatically  loads  the 
digital  state  feedback  control  subroutine  SFBCON  (see  Appendix  K). 
Rather  than  use  an  actual  microprocessor  and  associated  interface 
equipment,  it  was  decided  that  a  small  portion  of  the  PDP-15's 
memory  would  be  used  to  "simulate"  a  microprocessor.  Considering 
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the  present  state  of  microprocessor  technology,  this  substitution 
can  be  made  without  loss  of  applicability.  Following  the  initial 
setup,  program  control  is  transferred  to  the  logic  and  switching 
section  of  the  analog  computer.  SFBCON  runs  constantly  in  an 
endless  loop,  sampling  the  state  variables  at  a  rate  of  about 
71  samples  per  second  and  calculating  feedback  control  signals 
based  on  each  set  of  sampled  data.  Periphery  devices  include 
a  functional  miniature  joystick  and  an  analog  generated,  pictorial 
output,  both  used  for  pilot-in-the-loop  studies.  Of  course, 
in  addition  to  the  real-time  mode,  the  system  response  may  be 
observed  using  the  hybrid's  high  speed  repetitive  operation 
(rep-op)  mode.  Push-button  logic  allows  instantaneous  switching 
between  the  four  feedback  options:  (1)  no  feedback,  no  SAS; 

(2)  present  SAS;  (3)  continuous  (analog)  state  feedback;  and 
(4)  digital  state  feedback.  User  documentation  and  system 
diagrams  are  contained  in  Appendix  I. 


CHAPTER  9 
SUMMARY  OF  RESULTS 

The  results  of  this  research  effort  can  be  described  in  terms 
of  the  project's  two  major  accomplishments:  The  successful  imple¬ 
mentation  of  a  proven  mathematical  model  of  the  AV-8B  on  both  a 
large  scale  digital  system  and  an  analog/hybrid  system,  and  the 
design  and  implementation  of  an  effective  SAS  for  that  model.  If 
it  can  be  conclusively  determined  that  the  present  simulation 
accurately  describes  the  dynamic  behavior  of  the  Harrier,  and  if 
it  can  be  shown  that  the  state  feedback  control  system  stabilizes 
that  model,  then  there  is  a  strong  logical  argument  that  this 
control  system  satisfies  the  goals  of  the  research.  This  summary 
of  results  will  focus  on  showing  that  the  newly  designed  digital 
state  variable  controller  fulfills  this  project's  stated  objective. 

The  analysis  of  results  must  begin  with  the  mathematical 
model,  for  if  it  is  invalid,  the  possibility  exists  that  any  results 
obtained  from  it  could  be  fruit  of  the  proverbial  poisoned  tree. 

Of  the  validity  of  the  final  AV-8B  model,  Lebacqx  [7,  p.9]  comments: 
"...it  was  judged  that  the  linear  model  represented  the  nonlinear 
time  histories  reasonably  well  at  these  flight  conditions."  It 
therefore  seems  reasonable  to  conclude  that,  providing  conditions 
remain  such  that  small  perturbation  theory  holds  (as  described  in 
Chapter  4),  the  Calspan  TM  No.  98  linearized  stability  and  control 
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derivatives  define  a  valid  mathematical  model  of  the  AV-8B. 

The  next  logical  concern  is  the  accuracy  with  which  the 
digital  and  hybrid  simulations  described  earlier  in  this  report 
represent  the  Calspan  model.  The  answer  to  this  question  is  illus¬ 
trated  by  comparison  of  the  graphical  outputs  shown  in  Part  I 
of  Appendix  L  with  the  Cal  span -genera ted  curves  presented  in 
Appendix  M.  The  plots  are  nearly  identical  in  every  respect. 

Also,  since  the  hybrid  model  has  been  shown  to  be  in  close  accoru 
with  the  NATS  digital  simulation,  it  can  be  concluded  that  bot. 
computer  simulations  describe  the  actual  behavior  of  the  Harrier 
in  flight. 

Having  reached  this  favorable  conclusion  concerning  the  fidelity 
of  the  current  computer  representations  of  the  unmodified  Harrier, 
attention  can  be  focused  on  the  performance  of  the  digital  state 
feedback  SAS.  Reference  to  Parts  II  through  IV  of  Appendix  L 
indicates  that  application  of  the  newly-designed  digital  controller 
results  in  stabilization  of  both  the  lateral  and  the  longitudinal 
modes  of  motion  of  the  model.  The  stabilization  is  evidenced  by 
the  fact  that  the  state  variables  are  driven  to  zero.  Recall  that, 
due  to  the  nature  of  the  linearization  process,  the  state  variables 
of  the  model  represent  disturbance  quantities.  Thus,  the  regulation 
of  these  variables  to  zero  has  the  physical  significance  of  returning 


the  aircraft  to  its  "trim"  condition,  and  therefore  of  stabilizing 
it.  A  closer  examination  of  the  time  responses  in  Appendix  L  is 
in  order. 

The  original  objectives  of  this  research  called  for  an 
analysis  of  the  intake  momentum  drag  instability  in  the  transition 
region.  By  setting  the  initial  conditions  of  the  model  to  reflect 
actual  values  the  state  variables  might  equal  at  the  outset  of  the 
" wea t here oc king"  phenomena.  Parts  II  through  IV  of  Appendix  L 
provide  an  interesting  comparative  study  of  the  performance  of  the 
present  SAS  versus  that  of  the  state  feedback  controller  in  this 
situation.  Parts  II  and  III  demonstrate  the  system's  effectiveness 
at  various  velocities  in  the  transition  region.  Part  II  treats 
the  lateral  modes,  while  Part  III  contains  longitudinal  responses. 

The  curves  of  Part  IV  present  a  progression  of  response  characteristics, 
from  that  of  the  unstable  bare  airframe  to  the  completely  stable, 
near  first  order  response  of  the  state  feedback  equipped  system. 

As  a  further  demonstration  of  the  ability  of  the  system  to 
recover  the  aircraft.  Parts  V  and  VI  of  Appendix  L  illustrate  the 
expansion  of  the  flight  envelope  that  implementation  of  the  digital 
SAS  allows.  As  discussed  earlier,  the  flight  manual  prohibits 
angle  of  attack  exceeding  15  units  and  side  velocities  exceeding 
30  kts.  The  state  feedback  controller  permits  these  limits  to 
be  safely  expanded,  at  least  as  far  as  the  linearized  model  is 
concerned. 
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A  thorough  study  of  the  time  trajectories  of  Appendix  L 
allows  one  to  formulate  several  generalizations  regarding  the 
relative  performance  of  the  two  stability  augmentation  packages. 

The  lateral  motion  of  the  aircraft  (presented  primarily  in  Part  II 
of  Appendix  L)  appears  to  be  close  to  neutrally  stable.  In  addition, 
the  relative  stability  of  the  aircraft  seems  to  increase  with 
velocity,  a  finding  not  in  keeping  with  the  predictions  of  earlier 
chapters  that  the  aircraft's  motion  was  most  unstable  between  30 
and  90  KIAS.  Despite  this  discrepancy,  ,it  was  judged  that  the  curves 
were  a  fairly  accurate  description  of  the  Harrier  with  existing 
stability  augmentation.  The  trajectories  of  the  state  feedback 
equipped  aircraft,  in  contrast  to  those  just  studied,  show  a  stable, 
nearly  first  order  response  in  all  cases.  Furthermore,  the  gain 
scheduling  algorithm  is  seen  to  do  an  adequate  job  in  maintaining 
a  damping  ratio  greater  than  0.5.  Recall  that  since  the  design 
locations  for  the  system  poles  lie  on  the  real  axis,  the  desired 
damping  ratio  is  unity.  Also,  note  how  intolerable  the  system  is 
to  large  disturbances  --  it  exhibits  a  "willingness"  to  reach 
comparatively  large  instantaneous  values  of  angular  velocities 

i 

in  order  to  eliminate  a  large  side  velocity,  for  example.  (See 
Parts  II  and  VI  of  Appendix  L)  This  tenacity  to  remain  at  or  near 
trim  is  a  convenient  characteristic  as  far  as  the  model  is  concerned, 
for  it  also  implies  that  the  system  seeks  to  remain  in  the  linear 
region  (where  the  linearized  derivatives  are  valid). 


t 
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The  longitudinal  response  characteristics  (highlighted  in 
Part  III  of  Appendix  L)  of  the  model  in  many  respects  resemble  those 
of  the  lateral  case.  Under  conventional  stability  augmentation, 
the  aircraft  experiences  a  very  slight  divergence  in  pitch  at  all 
velocities  in  the  VSTOL  regime.  Also,  a  pitch  distrubance  causes 
a  rapid  reduction  in  forward  velocity  (U)  in  all  cases.  As  in  the 
lateral  case,  the  motion  of  the  state  feedback  equipped  system  is 
characterized  by  stable  responses  of  high  damping  ratio. 

Preliminary  pilot-in-the-loop  handling  qualities  studies 
conducted  on  the  hybrid  computer  demonstrate  the  superiority  of  the 
state  feedback  system.  The  pilot  evaluation  task  was  defined  as 
maintaining  steady  transition-regime  (usually  30  KIAS)  flight 
down  a  100'  x  100'  corridor.  The  pilot  controlled  the  aileron  and 
elevator  via  a  miniature  joystick.  His  performance  was  observed 
on  an  oscilloscope  using  an  analog-generated  pictoral  output.  The 
experimenter  had  the  ability  to  introduce  wind  gusts  of  varying 
magnitudes  from  any  direction.  The  pilot  was  seldom  able  to  maintain 
control  of  the  aircraft  under  conventional  stability  augmentation 
for  more  than  several  seconds.  On  the  other  hand,  the  state 
feedback  system  allowed  the  pilot  to  retain  control  for  an  extended 
period  of  time,  even  in  the  presence  of  turbulence.  The  pictoral 
responses  of  both  the  handling  qualities  studies  and  of  the  open  loop 
studies  were  video-taped  to  facilitate  further  study. 


As  stated  earlier,  the  objective  of  this  research  has  been 
to  design  an  improved,  microprocessor-compatible  three-axis  SAS 
for  the  AV-8B  Harrier.  Further,  the  SAS  was  to  be  capable  of 
functioning  at  two  levels  of  control  authority  -  a  "normal"  mode 
and  a  "recovery"  mode.  Closed  loop  pole  placement  by  digital 
state  feedback  is  a  suitable  method  to  meet  these  criteria.  The 
performance  of  the  "recovery"  mode  is  demonstrated  in  this  paper. 

Note  again  the  tenacity  with  which  the  system  seeks  to  return  the 
aircraft  to  a  trimmed  condition.  Further,  this  level  of  performance 
is  achieved  without  exceeding  the  control  authority  limit  imposed 
by  military  specifications.  [10]  There  are  two  possible  avenues 
that  may  be  followed  in  implementing  a  "normal"  control  authority 
mode:  either  to  define  new  design  pole  locations  and  recalculate 
the  feedback  gain  matrices,  or  to  simply  reduce  the  control  authority 
of  the  present  SAS.  Preliminary  studies  indicate  that,  from  a 
handling  qualities  point  of  view,  the  latter  approach  is  preferable. 


CHAPTER  10 


CONCLUSIONS 


The  general  conclusion  that  may  be  drawn  from  this  study  is  that 
it  is  feasible  to  stabilize  the  motion  of  the  AV-8B  Harrier  in  low 
speed  flight  by  using  a  microprocessor-based  state  feedback  controller 
to  "directly"  relocate  the  system  poles.  This  controller  design 
has  been  shown  to  be  capable  of  fully  stabilizing  a  proven  computer 
simulation  of  the  Harrier.  The  achievement  of  these  objectives 
has  required  the  accomplishment  of  the  following  tasks: 


•  The  computer  implementation  and  verification 
of  the  Cal  span  linearized  mathematical  model 
for  the  AV-8B  Harrier  airframe.  The  computer 
model  includes  a  least-squares  interpolation 
routine  to  accommodate  the  assignment  of  the 
model's  variable  stability  gains.  The  research 
required  the  development  of  both  a  digital 
and  a  hybrid  computer  simulation.  The  models 
were  verified  by  comparison  of  time  trajectories 
to  those  in  the  literature. 


•  The  construction  and  computer  implementation  of 
a  valid  model  for  the  existing  Harrier  SAS.  The 
model  required  the  use  of  information  from  various 
sources  in  the  literature. 


t  The  development  of  a  computer  algorithm  for  the 
direct  relocation  of  the  closed  loop  system  poles. 
The  program  was  then  used  to  determine  the  lateral 
and  longitudinal  state  feedback  matrices  required 
to  meet  the  desired  pole  location  criteria  over  the 
studied  velocity  range  (0-100  KIAS). 


•  The  design  of  a  microprocessor-compatible  controller 
to  implement  and  test  the  calculated  state  feedback 
gains. 
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•  A  comparison  of  the  performance  of  the  existing 
conventional  SAS  design  versus  that  of  the  digital 
device. 


A  number  of  conclusions  concerning  the  simulation  of  VSTOL 
aircraft,  digital  flight  control,  and  the  relative  performance  of 
the  two  Harrier  stability  augmentation  packages  may  be  formulated 
based  on  this  research.  The  statements  made  regarding  flying 
qualities  are  the  result  of  preliminary,  informally  conducted  studies. 
The  conclusions  are  summarized  below,  and  additional  details  are 
contained  in  the  text. 


•  The  present  computer  simulations  represent  a 
fairly  accurate  and  versatile  tool  for  the 
design  and  evaluation  of  VSTOL  flight  control 
systems.  The  Naval  Academy’s  digital  and 
hybrid  computer  facilities  are  clearly  well 
adapted  to  the  study  of  such  systems. 


•  Direct  relocation  of  closed  loop  poles  using 
state  feedback  is  an  effective  method  of 
flight  control  design.  Consisting  of  a 
blend  of  modern  (state  feedback)  and  classical 
(root  locus)  control  theories,  this  procedure 
benefits  from  the  user's  ability  to  rapidly 
visualize  expected  results  --  a  feature  lacking 
in  many  control  design  methodologies. 

•  The  sampling  rate  allowed  by  this  controller  design 
(a  rapid  71  sps.)  is  more  than  adequate  to  insure 
performance  quite  comparable  to  that  of  an  analogous 
continuous  system. 


•  In  all  cases,  the  state  variable  controller  was 
found  to  stabilize  the  Harrier,  providing  clean, 
nearly  first  order  response.  It  is  judged,  based 
on  the  response  characteristics,  that  the  digital 
controller  facilitated  movement  of  the  system 
poles  to  locations  in  close  proximity  to  the  design 
values. 
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•  In  contrast,  the  existing  conventional  SAS, 
while  providing  definite  improvement  in  the 
aircraft's  handling  qualities  is  unable  to 
fully  stabilize  its  motion.  In  general, 
the  responses  of  the  conventional  SAS  model 
exhibit  a  near  neutral  stability. 


•  Pilot-in-the-loop  handling  qualities  studies 
demonstrate  clearly  the  superiority  of  the 
digital  state  feedback  system.  While  the 
pilot  is,  in  general,  unable  to  maintain 
control  of  the  conventionally-augmented 
system,  the  system  equipped  with  the  state 
feedback  device  is  quite  controllable,  even 
in  turbulent  conditions. 


•  The  lateral  acceleration  feedback  loop  of  the 
hybrid  simulation's  conventional  SAS  model 
proved  so  susceptible  to  the  effects  of 
spurious  noise  in  the  system  that  it  had  to  be 
removed  from  the  circuit.  This  result,  coupled 
with  similar  findings  found  in  the  literature, 
casts  some  degree  of  doubt  on  the  use  of  the 
actual  system  on  the  Harrier,  and  particularly 
on  its  very  light  detection  threshold  (.06G). 


The  determination  of  whether  the  controller  designed  in 
this  project  is  capable  of  stabilizing  the  Harrier  rests  on  a 
number  of  as  yet  unanswered  questions.  The  following  are  among 
the  key  considerations: 


•  The  effect  of  many  of  the  higher  order  terms 
(such  as  instrumentation  delays,  more  complete 
actuator  dynamics  models,  etc.)  that  have  been 
omitted  from  this  model. 


•  The  limitations  associated  with  the  use  of  a 
linearized  model.  As  shown  in  Chapter  4, 
linearization  of  the  stability  derivatives 
decouples  the  lateral  and  longitudinal  equations 
of  motion.  However,  our  discussion  of  the 
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•  intake  momentum  drag  phenomena  alluded  to 

probable  coupling  between  the  lateral  and  longi¬ 
tudinal  modes.  The  strength  of  this  coupling 
and  its  effect  on  the  controller  have  yet  to  be 
determined. 


•  The  effects  of  sampling  rate,  pure  time  delays, 
and  quantization  errors.  A  discussion  of  on¬ 
going  research  in  regard  to  these  and  other 
important  parameters  in  the  study  of  digital 
flight  control  is  contained  in  reference  16. 


•  The  ability  to  instrument  the  aircraft  to  measure 
all  necessary  system  states.  The  use  of  an 
observer  may  be  required  to  approximate  unavailable- 
system  states. 


•  The  controller's  ability  to  cope  with  model 
inaccuracies  and  to  adapt  to  inevitable  changes 
in  plant  dynamics.  Complete  sensitivity  and 
robustness  studies  are  clearly  in  order. 


As  the  study  of  the  state  feedback  control  of  aircraft  is  in  its 
infancy,  the  literature  provides  little  help  with  the  many  questions 
that  this  research  has  prompted. 

In  surmiary,  an  accurate  representation  of  the  dynamic  behavior 
of  the  AV-8B  Harrier  has  been  programmed  on  both  the  Naval  Academy 
PDP-15/EAI  681  Hybrid  Computer  System  and  on  the  Naval  Academy 
Time-Sharing  System  (NATS).  In  addition,  a  digitally-based  variable 
gain  state  feedback  controller  has  been  added  to  the  simulations, 
and  proved  to  be  capable  of  fully  stabilizing  them.  While  it  is 
doubtful  that  the  controller  in  its  present  form  is  the  final  answer, 
it  has  been  demonstrated  beyond  a  reasonable  doubt  that  microprocessor- 
based  variable  gain  state  feedback  control  is  a  valid  approach  to 
the  final  solution  of  the  Harrier's  stability  and  control  problems. 
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APPENDIX  A 
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THE  AV-8B  STABILITY  DERIVATIVES 

This  appendix  presents  the  stability  and  control  derivatives 
for  the  AV-8B  Harrier  linear  mathematical  model  formulated 
by  Cal  span  Corporation.  This  data  has  been  extracted  from 
TM  No.  98,  and  is  included  here  in  the  interest  of  completeness. 
The  A  and  B  matrices  shown  in  Figures  4.3  and  4.4  may  be 
derived  from  this  data  by  simple  substitution. 


AV-8B  LONGITUDINAL  MODEL 


STABILITY  AND  CONTROL  DERIVATIVES 


AV-8B  LATERAL  MODEL 


VQ( ft/sec) 


STABILITY  AND  CONTROL  DERIVATIVES 

0  50.668  84.447  109.781  135.116  177.339 


APPENDIX  B 

ROOTS  OF  THE  AV-8B  CHARACTERISTIC  EQUATION 

The  roots  of  the  system  characteristic  equation  are  tabulated 
and  plotted  as  a  function  of  velocity. 


ROOTS  OF  THE  AV-8B  LINEAR  AIRFRAME  MODEL 


LATERAL 


LONGITUDINAL 


AV-8B  ROOT  POSITION  AS  A 


APPENDIX  C 


MATRIX  MANIPULATION  SUBROUTINES 
( "MATPAK" ) 

Nearly  all  of  the  programs  utilized  in  conjunction  with  this 
research  rely  heavily  on  matrix  techniques.  The  file  listed  in 
this  appendix  contains  many  of  the  most  commonly  used  algorithms. 


I 


MATPAK 


05/05/79 


16:  51:  23 
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***  SUBROUTINE  TO  TAKE  THE  NEGATIVE  TRANSPOSE 

***  OF  A  MATRIX 

SUBROUTINE  NTRAN< X, K, L, Y, M,  N) 

DIMENSION  X<20. 20), Y<20, 20) 

DO  10  I  =  liK 
DO  10  J=1 »  L 
10  Y<J,  I)  =  -X(I,  J) 

M  =  L 
N  =  K 
RETURN 
END 

***  SUBROUTINE  TO  TAKE  THE  TRANSPOSE 

***  OF  A  MATRIX 

SUBROUTINE  TRAM  X,  K>  L,  Y,  M,  N) 

DIMENSION  X  <  20*  20 ) ,  Y  < 20,  20 ) 

DO  10  1=1, K 
DO  10  J= 1 , L 
10  Y(  J,  I )  =  X(  I,  J) 

M  =  L 
N  =  K 
RETURN 
END 

***  SUBROUTINE  TO  ADD  TWO  MATRICES 
SUBROUTINE  MADD< R,  A,  N,  M,  B) 

DIMENSION  R<  20, 20), A(20, 20),  B(20,  20) 

* 

DO  10  1=1, N 
DO  10  J= 1 , M 

10  R(  I,  U)  =  A(  I,  J)  +  B(  I,  J) 

RETURN 

END 

***  SUBROUTINE  TO  SUBTRACT  TWO  MATRICES 
SUBROUTINE  MSUB(R,  A,  N,  M,  F> 

DIMENSION  R<  20, 20), A<20, 20), B<20, 20) 

* 

DO  10  1  =  1,  N 
DO  10  J-l,  M 

10  R<  I,  J)  =  A(  I,  J)  -  B<  I,  J) 

RETURN 

END 

***  SUBROUTINE  TO  MULTIPLY  TWO  MATRICES 
SUBROUT  I  NEMMULT(R,  A,  N,  M,  L,  B) 

DIMENSION  R ( N,  L ) ,  A ( N,  M ) ,  B < M,  L ) 

DIMENSION  V<20,  20) 

* 

DO  10  1  =  1,  N 
DO  10  J=1,L 
V<  I,  J)  =  0. 

DO  10  K-l,  M 

V(I,J)  =  A<  I,  K)*B(K,  J)  +  V(I,J) 

10  CONTINUE 
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MATPAK 


DO  20  1=1, N 
DO  20  J=l,  L 
20  R( I,  J)  = 
RETURN 
END 


SUBROUTINE  MINV 


PURPOSE 

INVERT  A  MATRIX 


USAGE 


CALL  MINV<  A,  N,  D,  L,  M) 


DESCRIPTION  OF  PARAMETERS 

A  -  INPUT  MATRIX,  DESTROYED  IN  COMPUTATION  AND  REPLACED  BY 
RESULTANT  INVERSE. 

N  -  ORDER  OF  MATRIX  A 
D  -  RESULTANT  DETERMINANT 
L  -  WORK  VECTOR  OF  LENGTH  N 
M  -  WORK  VECTOR  OF  LENGTH  N 


REMARKS 

MATRIX  A  MUST  BE  A  GENERAL  MATRIX 


SUBROUTINES  AND  FUNCTION  SUBPROGRAMS  REQUIRED 
NONE 


METHOD 


THE  STANDARD  GAUSS- JORDAN  METHOD  IS  USED.  THE  DETERMINANT 
IS  ALSO  CALCULATED.  A  DETERMINANT  OF  ZERO  INDICATES  THAT 
THE  MATRIX  IS  SINGULAR. 


SUBROUTINE  MINV'R, X, N, D1 ) 

DOUBLE  PRECISION  A, D, B IGA,  HOLD 

DIMENSION  A< 100),  R< 10,  10 ) , L ( 100) ,  M ( 100) 
DIMENSION  X< 10, 10) 

K=0 

DO  9999  1  =  1,  N 
DO  9999  J=1,N 
K=K+1 

9999  A(K)=X< I, J) 


IF  A  DOUBLE  PRECISION  VERSION  OF  THIS  ROUTINE  IS  DESIRED,  THE 
C  IN  COLUMN  1  SHOULD  BE  REMOVED  FROM  THE  DOUBLE  PRECISION 
STATEMENT  WHICH  FOLLOWS. 


MATPAK 


05/05/79 


16: 52: 16 
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*  DOUBLE  PRECISION  A,  D,  B IGA,  HOLD 

*  SEE  ABOVE 

» 

*  THE  C  MUST  ALSO  BE  REMOVED  FROM  DOUBLE  PRECISION  STATEMENTS 

*  APPEARING  IN  OTHER  ROUTINES  USED  IN  CONJUNCTION  WITH  THIS 

*  ROUTINE. 

* 

*  THE  DOUBLE  PRECISION  VERSION  OF  THIS  SUBROUTINE  MUST  ALSO 

*  CONTAIN  DOUBLE  PRECISION  FORTRAN  FUNCTIONS.  ABS  IN  STATEMENT 

*  10  MUST  BE  CHANGED  TO  DABS. 


SEARCH  FOR  LARGEST  ELEMENT 

D=l.  0 
NK=-N 

DO  80  K=l, N 
NK=NK+N 
L  <  K ) =K 
M<K)=K 
KK=NK+K 
BIGA=A<  KK ) 

DO  20  J=K,  N 
IZ=N*( J-l ) 

DO  20  I=K, N 
IJ=IZ+I 

10  IF(  DABS(BIGA)-  DABS(A(IJ)))  15,20,20 
15  B I GA= A ( I J ) 

L<K)=I 
M(K)=J 
20  CONTINUE 

INTERCHANGE  ROWS 

J=L  <  K ) 

IF(J-K)  35,35,25 
25  KI=K-N 

DO  30  1=1, N 

KI=KI+N 

HOLD=-A(KI) 

JI=KI-K+J 
A<KI ) =A<  JI ) 

30  A( JI )  =HOLD 

INTERCHANGE  COLUMNS 

35  I =M  <  K ) 

IF(I-K)  45,45,38 
38  JP=N*  < I - 1 ) 


MK 


I 

t 
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DO  40  J=l, N 
JK=NK+J 
JI=JP+J 
HOLD=-A ( JK ) 

A< JK)=A< JI ) 

40  A<  JI )  =HOLD 

DIVIDE  COLUMN  BY  MINUS  PIVOT  (VALUE  :  /  V  El.  EMEr ' 

CONTAINED  IN  BIGA) 

45  IF  (BIGA)  48,46,43 

46  D=0.  0 

GO  TO  160 
48  DO  55  1=1, N 

IF(I-K)  50,55,50 
50  IK=NK+I 

A( IK)=A( IK) /(-BIGA) 

55  CONTINUE 

REDUCE  MATRIX 

DO  65  1  =  1,  N 
IK=NK+I 
HOLD=A  < IK) 

I J=I-N 
DO  65  J=l, N 
I J=I J+N 

IF(I-K)  60,65,60 
60  IF(J-K)  62,65,62 
62  KJ=I J-I+K 

A( IJ)=HOLD*A(KJ)+A( I J) 

65  CONTINUE 

DIVIDE  ROW  BY  PIVOT 

KJ=K— N 
DO  75  J=l,  N 
KJ=KJ+N 

IF(J-K)  70,75,70 
70  A(KJ)=A(KJ)/BIGA 
75  CONTINUE 

PRODUCT  OF  PIVOTS 

D=D*BIGA 

REPLACE  PIVOT  BY  RECIPROCAL 

A( KK )  =  1.  O/BIGA 
80  CONTINUE 


96 
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*  FINAL  ROW  AND  COLUMN  INTERCHANGE 

« 

K=N 

100  K=  <  K- 1 ) 

IF<K)  150,  150,  105 
105  I=L<K> 

IF(I-K)  120,120,108 
108  JQ=N*<K-1> 

JR=N*< 1-1 ) 

DO  110  J=l,  N 
JK=JQ+J 
HOLD=A  <  JK ) 

JI=JR+J 

A(JK)=-A(JI) 

110  A(JI)  =HOLD 
120  J=M<K> 

IF< J-K)  100,  100,  125 
125  KI=K-N 

DO  130  1=1,  N 
KI=KI+N 
HOLD=A(KI ) 

JI=KI-K+J 
A(KI )=-A< JI ) 

130  A( JI )  =HOLD 
GO  TO  100 
150  GO  TO  160 
160  K=0 
DO  170  1  =  1,  N 
DO  170  J=1,N 
K=K+1 

170  R(I,J)  =  A(K) 

D1=D 

RETURN 

END 


I  li 
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APPENDIX  D 


STATE  FEEDBACK  CONTROL  DESIGN  PROGRAM 


(" FADDM1 ") 


This  appendix  presents  FADDM1 ,  the  state  feedback  gains 
calculating  algorithm  and  its  accompanying  subroutines.  FADDM1 
was  the  major  controller  design  tool  of  the  research.  A  sample 
run  is  included  to  familiarize  the  reader  with  its  use,  and  to 
demonstrate  that  nr  unique  sets  of  feedback  gains  exist  for  each 
specification  of  root  locations,  where  n  is  the  system  order  and 
r  is  the  number  of  inputs. 


FADDM1 


05/05/79 


17:  37: 39 
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*  PROGRAMMER:  E.  E.  MITCHELL 

*  PROGRAM  NAME:  FADDM1 

*  SYSTEM:  DTSS  FORTRAN 

*  FUNCTION:  COMPUTES  STATE  FEEDBACK  GAINS  FOR 

*  THE  SOLUTION  OF  THE  ARBITRARY  POLE 

*  PLACEMENT  PROBLEM  FOR  MULT I -INPUT 

*  MULT I -OUTPUT  SYSTEMS 

* 

***  MAIN  PROGRAM  FOR  FADD 
* 

********************************************* 

* 

DIMENSION  A( 10, 10),  H( 10,  10). S( 10,  10) 

DIMENSION  P<  10),  AX<  10,  10),  B<  10,  10,  10) 

DIMENSION  C( 10,  10) 

DIMENSION  Q(10),RR<10),RI<10),D(10,  10,  10) 
CHARACTER  YES, FILE1.  FILE2 
LIBRARY "*£05006:  FADD" 

L I BRARY " *E05006 : MATPAK" 

L I BRARY " *E05006 :  PRNT" 

LI BRARY "*E05006:  EVAL" 

M  =  10 

IN=0 

I0UT=0 

*  DATA  FROM  FILE  OR  TTY 
PRINT,  "INPUT  FROM  A  FILE?" 

INPUT,  YES 

IF ( YES  .  NE.  "YES")  GO  TO  5 
PRINT,  "FILE  NAME  IS?" 

INPUT,  FI LEI 
OPENFILE  2, FILE1 
IN=2 


*  OUTPUT  TO  FILE  OR  TTY 

* 

5  CONTINUE 

PRINT, "OUTPUT  TO  A  FILE?" 

INPUT,  YES 

IF< YES  .  NE.  "YES")  GO  TO  10 
PRINT, "FILE  NAME  IS  ?" 

INPUT,  FILE2 
OPENFILE  3, FILE2 
REWIND  3 
ENDFILE  3 
101*3 

10  CONTINUE 

* 

*  READ  DATA  FROM  FILE  OR  TTY  DEPENDING  UPON  IN=0  OR  2 

*  IN*0  FOR  TTY,  IN*2  FOR  FI LEI 
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IF<  IN  EQ.  0 )  PR  I  NT 

IF<  IN  .  EQ.  0)  PR  I  NT/  "A  =  N  BY  N  MATRIX" 

IF(  IN  .  EQ.  0)  PR  I  NT 

IF  <  IN  .  EQ.  0)  PR  I  NT/  "N  =  " 

READ (IN,  500 )N 

IF ( IN  .  EQ.  0)PRINT,  "ENTER  A  BY  ROWS,  1  ROW/LINE" 

PRINT 

DO  8  1  =  1,  N 

8  READ< IN, 500) <A< I,  J), J=l, N) 

IF<  I  OUT  .  NE.  0 )  WR I TE  <  I  OUT  ,510) 

WR I TE  < 10, 510) 

DO  30  1  =  1,  N 

I F ( I OUT  .  NE.  0) WRITE < I OUT, 520) ( A(  I ,  J) ,  J=1 ,  N  > 

30  WRITE <  10,  520)  <A(I,J),J=1,N) 

CAUL  FADD( A,  B,  AX,  P,  N,  M,  101,  102) 

* 

« 

500  FORMAT < V) 

510  FORMAT <////“  A  MATRIX") 

520  FORMAT (  "0" »  5< G12.  5,  2X ) ,  5< /  "  ",  4(012  5,  2X  )  )  ) 

IF (  IN  .  EQ.  0)PRINT»  "N,  L,  B  BY  ROWS" 

READ ( IN,  500) N,  L 
DO  40  1  =  1,  N 

40  READ  (IN,  500  )(C(I,J),J=1,L) 

DO  45  1=1, N 
DO  45  J=1 , L 
45  H ( I »  J ) =C ( I ,  J ) 

CALL  PRM(H,  "INPUT  B",  10,  10,  N,  L,  0) 

DO  70  K=l, N 
DO  50  1=1, N 
DO  50  J=l, N 
50  AX ( I, J)=B(K,  I,  J) 

CALL  MMULT(AX,  AX,  N,  N,  L,  C) 

DO  60  1=1, N 
DO  60  J= 1 , L 
60  B(K,  I,  J)=AX< 1,  J) 

70  CONTINUE 

IF( IN  .  EQ.  0)PRINT, "ENTER  N  ROOT  LOCATIONS,  RR , R 1 ,  1/LINE 
DO  80  1  =  1,  N 

80  READ (IN, 500 ) RR ( I ) , R I ( I ) 

PRINT, "DESIRED  ROOT  LOCATIONS" 

DO  85  1  =  1,  N 

85  PRINT  520, RR(I),RI(I) 

DO  100  Ll=l , N 
*  CALL  PRRV(P, "P",  N+l,  0) 


PR1=RR(L1 ) 
PI 1=RI (LI ) 
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90  Q<K)=B<K,  I.  J) 

*  CALL  PRRV<Q,  "Q",  N,  0) 

PR=RR<L1 ) 

PI=RI<L1) 

CALL  EVAL<Q,  N,  PR,  PI ) 

D<  LI,  I, J ) =PR/PR 1 

*  PRINT,  D<L1,  I,  J),  PR,  PR1 
100  CONTINUE 

DO  110  1*1,  L 
DO  110  J=l,  N 
C(I,  J )  =0. 

110  C(l,  J)-l. 

DO  180  Kk=l,L*U 
DO  130  1=1,  N 
DO  130  J=1,L 
S(  J.  I  )=-C(  J,  I ) 

IF<C(J,  I)  .  EQ.  0.  )GQ  TO  130 

DO  120  K=1,N 

120  AX  ( K,  I )  =  D<  I,  K,  J) 

130  CONTINUE 
1  =  1 

150  DO  160  J=l» L 

IF(C(J,  I)  .  EQ.  0)G0  TO  180 

C<J,  I )  =0. 

IF .  EQ.  L)GO  TO  170 
C<J+1,  I)  =  l. 

GO  TO  175 
160  CONTINUE 
170  C<1,  I)=l. 

1  =  1  +  1 

IF<  I  .  GT.  N)GO  TO  175 
GO  TO  150 
175  CONTINUE 

*  CALL  PRMtC,  "C",  10,  10,  N,  N,  0) 

*  CALL  PRM(AX,  "AX",  10,  10,  N,  N,  0) 
CALL  MINV<AX,  AX,  N,  DET) 

IF<DET  .  NE.  0)  GO  TO  190 
PRINT, "DET  =  0" 

CALL  PRM(S,  "C",  10,  10,  N,  L,  0) 

GO  TO  180 

190  CALL  MMULT(AX,  S,  L,  N,  N,  AX) 

CALL  PRM<  AX,  "K  GAINS",  10,  10,  L,  N,  0) 
CALL  MMULT(AX,H,  N,  L,  N,  AX) 

CALL  PRM<  AX,  "B*K ",  10,  10,  N,  N,  0) 
CALL  MSUB<AX,  A,  N,  N,  AX) 

CALL  PRM(  AX,  "A-BK",  10,  10,  N,  N,  0) 

IF(  103.  NE.  3)  GO  TO  180 
J  CALL  PRM<  AX,  "A-BK",  10,  10,  N,  N,  103) 

180  CONTINUE 
.  END 
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SUBROUTINE  FADD(  A,  B,  AX,  P,  N,  M,  10,  102) 


*  SUBROUTINE  TO  FIND  THE  INVERSE  OF  THE  (SI -A) 

*  MATRIX  -  A  MATRIX  IN  GENERAL  FORM 


E.  E.  MITCHELL 


*  CALL  VARIABLES 


6/78 


*  A  =  A  MATRIX  IN  (SI-A) 

*  B  =  COEFFICIENTS  OF  NUMERATOR  POL YNOM INALS 

*  IN  THE  (SI-A)  INVERSE 

*  B(  1,  I,  J)*S**(N-1 )  +  B(2,  I,  J)*S**(N-2)  +. 

*  B(N,  I,  J) 

*  P  =  VECTOR  OF  CHARACTERISTIC  EQN  COEFFICIENTS 

**  P(1)*S**N  +  P(2)*S**(N-1 )  +.  .  .  +  P(N) 

*  AX  =  INVERSE  OF  A  MATRIX  UPON  RETURN 

*  B ( N+ 1 ,  1 ,  1 )  =  DETERMINANT  OF  A  MATRIX  UPON  RE'-RN 

*  M  =  OR  >  N+l ,  EQUALS  THE  DIMENSION  OF  A, B, AX  &  F 

*  MATRICES  IN  THE  CALLING  PROGRAM 

*  N  =  ORDER  OF  A  MATRIX 

*  101  =  OUTPUT  FILE  #  FOR  INVERSE  OF  (SI-A) 

*  102  =  OUTPUT  FILE  #  FOR  INVERSE  OF  A  MATRIX 


*#*#*##**##*#*******#*****##<HHHHHHHHHt#1t#1Hm-jHt 

DIMENSION  P(M) ,  A(M,  M) ,  AX  < M,  M) ,  B< M,  M,  M ) 

L I BRARY " *E05006 :  MATPAK" 

15  DO  20  1  =  1,  N 
DO  20  J=l, N 
20  B(  1,  I,  J)=0.  0 
DO  25  1=1, N 
25  B(l,  I,  I)  =  1.  0 
* 

*  MAKE  A1  =  A 

* 

DO  35  1=1, N 
DO  35  J=l, N 
35  AX ( I, J)  =  A(I, J) 

P(l>  =  -1.  0 

* 

*  NEED  A  COUNTER,  K=2, N+l  FOR  THE  N  MATf  ICES 


*  THE  LOOP  AROUND  THE  ALGORITHM 

* 

50  K=K+1 
XK=K— 1 
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*  FIRST  POO 

*  PRINT, "Gl" 

* 

POO  »  0. 

DO  60  1*1,  N 

60  POO  =  POO  +  AX  ( I,  I) 

T0L=1.  E-6 
P<K)  *  POO/XK 

*  PRINT, " P ( K ) " ,  K,  P(K) 

*  THE  B  MATRICES 
DO  80  1*1,  N 

DO  80  J*l,  N 

B(K,  I,  J)  =  AX  <  I,  <J)  -  P<K)*B< 1,  I,  J) 

80  CONTINUE 

*  WRITE < 10,  620 )K 

*  DO  82  1*1,  N 

*  82  WRITE<  10,  560)  (BOC,  I,  J),  J*l,  N) 

DO  86  1*1,  N 

DO  86  J=1 , N 

86  IF(ABS<B(K,  I,  J)  )  .  LE.  TOL)  B<K,  I,  J)  =0.0 

*  THEN  THE  A I  MATRICES 
IF<K  .  EQ.  N+l )  GO  TO  120 
DO  100  1=1,  N 

DO  100  J=l»  N 
AX  ( I,  J)  =  0.  0 
DO  100  L=1 ,  N 

AX  ( I,  J)  =  AX(I.J)  +  A<  I,  L)*B<K,  L,  J) 

100  CONTINUE 

*  WRITE <10,  622) 

*  DO  103  1  =  1,  N 

*  103  WRITE< 10,  560) <AX< I,  J),  J*l,  N) 

DO  105  1  =  1,  N 

DO  105  J=1 , N 

105  IF<ABS<AX<  I,  J)  )  .  LE.  TOL)  AX<I,J)  =0.0 
GO  TO  50 

*  END  OF  LOOP 
120  WRITE<IO,  505) 

ISW=0 

DO  140  1*1,  N 

WRITE< 10,  560) < B < N+ 1 ,  I, J ) >  J= 1 ,  N) 

140  IF<B<N+1,  I,  J)  .  GT.  1.  E-05)ISW  =  1 
I F <  I SW  .  EQ.  0)  GO  TO  160 
PRINT, "ST0P=1,  CONTINUED" 

IF<  IN  .  EQ.  0)G0  TO  160 
INPUT,  ISW 

IF <  ISW  .  EQ.  1 )  STOP 
160  CONTINUE 

* 

*  HERE  TO  OUTPUT  RESOLVANT  MATRIX 

« 

WRITE< 10, 510)N 
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DO  180  1=1,  N+l 
J  =  N+l-I 
P!I)=-P!I) 

180  WRITEdO,  520)P!I),  J 
WRITE! 10,  530) 

DO  200  K=l,  N 
WRITE! 10,  550  >K 
DO  200  1  =  1,  N 

200  WRITE! 10/560) !B!K,  I,  J),  J=l,  N> 

GO  TO  285 
210  CONTINUE 
WRITE! 10,  570) 

WRITE! 10,  580) 

WRITE! 10.  585 )N,  N 
DO  220  1  =  1,  N 
DO  220  J=1,N 
WRITE! 10,  580) 

WRITE! 10,  590)1,  J 
DO  220  K=2, N+l 
L  =  N+l— K 

220  WRITE! 10,  595)B!K,  I, J), L 
250  CONTINUE 

*  OUTPUT  THE  INVERSE  OF  A 
260  XK  =— ! ( — 1 ) ** ! N— 1 ) ) *P ! N+ 1 ) 

B!N+1,  1,1)  =  XK 
WRITEdO,  60OXK 
WRITE! 10,  610) 

XK=-P!N+1 ) 

DO  270  1  =  1,  N 

WRITE! 10,  560) !B!N,  I,  J)/XK,  J=l,  N) 

DO  280  J= 1 , N 

280  AX ! I , J )  =  -B!N,  I,  J)/P!N+1 ) 

270  CONTINUE 
285  RETURN 
#* 

** 

** 

500  FORMAT!"  ",  6!G8.  2,  2X  )/"  ",  3X,  5!G8.  2,  2X  )  ) 

505  FORMAT!"!",////,  "O  SOLUTION  CHECK,  IS  B=0  ?•' ) 
510  FORMAT!///"  CHARACTERISTIC  EQN,  ORDER=", 14) 
515  FORMAT!///"  INPUT  1  FOR  RESOLVENT  MATRIX" 
it  /"  INPUT  2  FOR  INVERSE  OF  A  MATRIX" 

/"  INPUT  0  FOR  BOTH") 

520  FORMAT!"  ",  G14.  5,  "*S**",  12) 

530  FORMAT!///"  ADJOINT  MATRICES,  B!I),  WHERE" 

?<  /"B!0)*S**N  +  B!  1  )#S**!N-1 )  +.  .  .  +  B!N-1)*S" 
it  "  +  B !  N )  "  ) 

550  FORMAT!///"  B!",I2, ")",) 

560  FORMAT !  "0" ,  5! G12.  5,  2X  ) ,  5! /  "  ",  4(G12.  5,  2X  )  )  ) 

570  FORMAT!///"  RESOLVENT  MART I X  POL YNOM INALS" ) 
580  FORMAT!"  0") 
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585  FORMAT  <  "  +  I(S)  IS  A  ",  12,  "  BY",  12,  " 
590  FORMAT <  "+  /<  ",  12,  ",  ",  12,  "  )  ",  ) 

595  FORMAT  ("0",  4(011.  4,  "*S**“,  12,  2X  ) , 

&  5 ( / "  ",  3<G11.  4,  "*S**",  12,  2X)  )  ) 

600  FORMAT (///"  DETERMINANT  OF  A  =  ", G14. 
610  FORMAT ("  A  INVERSE") 

620  FORMAT ( "  B, K", 12) 

622  FORMAT <"  A,  K",  12) 

* 

END 
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SUBROUT  I NE  EVAL  !  B,  N,  Q,  S> 

* 

*  PROGRAMMER:  E.  E.  MITCHELL 

*  PROGRAM  NAME:  EVAL 

*  SYSTEM:  DTSS  FORTRAN 

*  FUNCTION:  POLYNOMIAL  EVALUATION  AT 

*  A  COMPLEX  POINT: 

*  A! 1 )+A!2)*X+A!3)*X**2+.  .  .  A!N+1  )*X-**N 

* 

*********#*********•****•»**#**•»••»**********#■«■•»■*«•#* 

* 

DIMENSION  A(10) 

DIMENSION  B< 10) 

COMPLEX  X,  R 

*  PRINT, “ENTER  POINT,  REAL  PART,  IMAG  PART" 

*  INPUT,  Q,  S 
DO  5  1=1, N 

5  A  < I ) =B  <  N+ 1  —  I ) 

X=CMPLX(Q,  S) 

R=CMPLX  ! 0.  ,  0.  ) 

J=N 

1  IF <  J)3,  3,  2 

2  R=R*X+A!J) 

J=.J- 1 

GO  TO  1 

3  CONTINUE 

*  PRINT  500, Q,  S 

*  PRINT  510,  R 
Q=REAL!R> 

S=A I MAG ( R ) 

RETURN 

12  CONTINUE 
G=SQRT ( Q*Q+S*S ) 

IF!  Q  .  EQ.  0.  0  .  AND.  S  .  EQ.  0.  0)G0  TO  10 
S=ATAN2  <  S,  REAL  !  R  >  >  *1 80.  /3.  141 5926 
PRINT  520,  Q,  S 
10  CONTINUE 
RETURN 

* 

500  FORMAT! “OVALUE  OF  POL YNOM INAL  AT  ",2G14.  5) 

510  FORMAT!”  RR,  RI  ",  2G14.  5) 

520  FORMAT !  "0  OR  MAGNITUDE  =",G14  5, /4X,  "ANGLE  G14.  5) 

END 
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RUN 

INPUT  FROM  A  FILE??  YES 
FILE  NAME  IS??  LON5 
OUTPUT  TO  A  FILE??  YES 
FILE  NAME  IS  ??  F 


A  MATRIX 


-.44000E-01 

-.57910E-03 

-32.160 

1.3320 

.72470E-03 

- .  381 60E-01 

-1 .6900 

-9.6633 

.OOOOOE-38 

.OOOOOE-38 

•OOOOOE-38 

1 .0000 

13680E-02 

.  459 1 0E-02 

.OOOOOE-38 

- . 66380E-01 

-0 .15940 
-0.34610 
.OOOOOE-38 
0.23080 
DESIRED  ROOT 

INPUT  D  MATRIX 
0.3371  1 
-2.5380 
.OOOOOE-38 
-. 36000E-01 
LOCATIONS 

-2.0000 

.OOOOOE-38 

-2.2000 

.OOOOOE-38 

-2.4000 

.OOOOOE-38 

-3.0000 

.OOOOOE-38 

-239.15 

.OOOOOE-38 

K  CAINS  MATRIX 
5396.0 
.OOOOOE-38 

52881 . 
.OOOOOE-38 

7967.4 

.OOOOOE-38 

38.120 
82.769 
. OOOOOE-38 
-55.195 

B*K  MATRIX 
-860.13 
-1867.6 
.OOOOOE-38 
1245.4 

-8429-3 

-18302. 

.OOOOOE-38 

12205. 

-1270.0 
-2757.5 
.  OOOOOF.-38 
1838.9 

-38.164 

-82.768 

.OOOOOE-38 

55.194 

A-BK  MATRIX 
860 .12 

1867.5 

.OOOOOE-38 

-1245.4 

8397.1 

18301 . 

. OOOOOE-38 
-12205- 

1271.3 

2747.9 

1 .0000 
-1839.0 
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K 

GAINS  MATRIX 

-2.0055 

-0 . 34484 

80.280 

30.519 

•33202E-01 

-0.75553 

-7.3926 

-1.1112 

B*K  MATRIX 

0.33086 

-0.19973 

-15.289 

-5.2393 

0.60983 

2.0369 

-9.0225 

-7.7422 

.OOOOOE-38 

•OOOOOE-38 

.  OOOOOE-38 

. 00000E 

-0.46406 

- . 52391 E-01 

18.795 

7.0837 

A-BK  MATRIX 

-0.37486 

0 .19975 

-16.871 

6.5713 

-0.609 10 

-2.0751 

7.3325 

-1 .921  1 

.OOOOOE-38 

.OOOOOE-38 

.OOOOOE-38 

1 .0000 

0.46269 

.56982E-01 

-18.795 

-7.1501 

K 

GAINS  MATRIX 

-1.8190 

-0 . 33227 

75.516 

29.800 

•36604E-01 

-0.83229 

-8.1443 

-1 .2242 

B*K  MATRIX 

0.30228 

-0.22761 

-14.783 

-5.1628 

0.53664 

2.2273 

-5.4660 

-7.2066 

.OOOOOE-38 

.OOOOOE-38 

.OOOOOE-38 

. OOOOOE 

-Q  .  4  2 1 1 3 

-.46727E-01 

17.722 

6.9218 

A-BK  MATRIX 

-0.34628 

0.22703 

-17.377 

6.4948 

-0.53592 

-2.2655 

3.7760 

-2.4567 

.OOOOOE-38 

.OOOOOE-38 

.OOOOOE-38 

1  .0000 

0.41977 

.5131 8  E— 0 1 

-17.722 

-6.9882 

K 

GAINS  MATRIX 

-4.7658 

-0.69854 

149-33 

40.941 

-20.043 

-3-3280 

494.84 

74.695 

B*K  MATRIX 

-5.9970 

-1 .0106 

143.01 

18.654 

52.518 

8.6882 

-1307.6 

-203.75 

.OOOOOE-38 

.OOOOOE-38 

.  OOOOOE-38 

.OOOOOE 

-0 .378*10 

- .  4 1 4  1 5  E  -  0 1 

16.652 

6.7603 

A-BK  MATRIX 

5.9530 

t.0100 

-175.17 

-17-322 

-52.5T8 

-8.7264 

1305.9 

194.08 

.OOOOOE-38 

.OOOOOE-38 

•OOOOOE-38 

1 .0000 

0.37703 

.46006E-01 

-16.652 

-6.8266 
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GAINS  MATRIX 

-1 .6640 

-0.32584 

71  .081 

29-050 

. 40005E-01 

-0.90904 

-8.8961 

-1.3373 

B*K  MATRIX 

0.27872 

-0.25451 

-14.329 

-5.0813 

0.47436 

2.4199 

-2.0229 

-6.6602 

.OOOOOE-38 

. OOOOOE-38 

•OOOOOE-38 

. 00000E 

-0.38548 

- .  42478E-01 

16.726 

6.7528 

A-BK  MATRIX 

-0.32272 

0.25393 

-17.831 

6.4133 

-0 . 4736*1 

-2.4581 

0.33289 

-3.0031 

.OOOOOE-38 

.OOOOOE-38 

.OOOOOE-38 

1 .0000 

0.38411 

•47069E-01 

-16.726 

-6.8192 

K 

GAINS  MATRIX 

-6.7910 

-0.94612 

206.31 

50.673 

-33-072 

-4.9150 

864.41 

138.31 

B*K  MATRIX 

-10.066 

-1 .5061 

258.52 

38.548 

86.286 

12.802 

-2265.3 

-368.57 

. OOOOOE-38 

.OOOOOE-38 

•OOOOOE-38 

. OOOOOE 

-0.37678 

-.41425E-01 

16.496 

6.7161 

A-BK  MATRIX 

10.022 

1.5055 

-290.68 

-37.216 

-86.285 

-12.840 

2263.6 

358.90 

•OOOOOE-38 

•OOOOOE-38 

.OOOOOE-38 

1 .0000 

0.37541 

.46016E-01 

-16.496 

-6.7825 

K 

GAINS  MATRIX 

-17-114 

-2.1710 

502.91 

102.28 

-97.707 

-12.583 

2723.1 

461  .99 

B*K  MATRIX 

-30.210 

-3.8956 

837.83 

139.44 

253-90 

32.686 

-7085.3 

-1207.9 

. OOOOOE-38 

.OOOOOE-38 

.  OOOOOE-38 

.OOOOOE 

-0.43246 

-. 48088E-01 

18.039 

6.9755 

A-BK  MATRIX 

30.166 

3-8951 

-869.99 

-138.10 

-253-90 

-32.724 

7083-6 

1198.3 

.OOOOOE-38 

. OOOOOE-38 

.OOOOOE-38 

1 .0000 

0.43109 

. 52679E-01 

-18.039 

-7.0419 

-109- 


-0.11082 

7.7953 


2.6455 
-19.746 
. 00000E 
-0.30621 


-2.6895 

19.747 

.OOOOOE 

0.30484 


K  GAINS  MATRIX 

-0.17151  -3-8712  11.223 

-0.17620  -421.37  -103-03 


B»K  MATRIX 
- . 3206 1 E-01  -141.43 

0.50657  1070.8 

38  . OOOOOE-38  . OOOOOE-38 

-.33242E-01  14.276 


-36.521 

257.61 

.OOOOOE-38 

6.2994 


38 


A-BK  MATRIX 


•31482E-01 

-0.54473 

.OOOOOE-38 

•37833E-01 


109.27 
-1072.5 
. OOOOOE-38 
-14.276 


37.853 
-267 .27 
1 .0000 
-6.3657 
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APPENDIX  E 

CALCULATION  OF  AV-8B  STABILITY  DERIVATIVES 
("MATCOF") 

The  subroutine  contained  in  this  appendix  calculates  the 
AV-8B  Harrier  stability  and  control  derivatives  at  any  velocity, 
V.  The  derivatives  are  returned  in  the  arrays  ALON,  BLON,  ALAT, 
and  BLAT.  The  subroutine  first  reads  in  the  coefficients  for  the 
predetermined  least  squares  polynomials.  Each  stability 
derivative  is  the  solution  of  a  polynomial  in  V,  the  input 
velocity. 
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SUBROUTINE  MATCOF < V,  ALON.  BLON,  ALAT.  BLAT) 

* 

♦PROGRAMMER:  R.  V.  WALTERS 
♦FILE  NAME:  MATCOF 
♦SYSTEM:  DTSS  FORTRAN 

♦FUNCTION:  SUBROUTINE  FOR  CALCULATION  OF 

♦  .  COEFFftS  OF  A  AND  B  MATRICES  FOR 

♦  AV-8  HARRIER  SIMULATION  AT  ANY  V. 

««*****»«***♦**•****»***#« 

♦ 

DIMENSION  AAL0NU6,  6),  BBL0N(9,  4 ) ,  AALAT (11,4),  BBLAT < 6,  4> 
DIMENSION  V5<6,  1 ) ,  V3<4, 1 ) 

DIMENSION  ALQN<16,  l ) ,  ALAT  <  1.6;  1>,BL0N<12,  1>,BLAT<12,  1) 

* 

LIBRARY  "MATPAK" 

♦ 

V5<  6,  1  )=1. 

DO  1  1=5*  1,  -1 

V5< I*  1 >=V5< 1*1,  D*V 

1  CONTINUE 
♦ 

V3(4,  1)=1. 

DO  2  1=3,  1,-1 

V3< I  *  1 >  =V3 (1*1,  1)^V 

2  CONTINUE 
♦ 

OPENFILE  6.  "ALON" 

REWIND  6 
* 

DO  4  J=l,  6 

DO  3  1  =  1,  16 

READ  <6,  >AALON< I,  J) 

5  CONTINUE 
4  CONTINUE 
» 

OPENFILE  7,  "BLON" 

REWIND  7 
♦ 

DO  7  J=l,  4 
DO  8  1=1,9 
READ  <7,  )BBLON(  I,  J) 

8  CONTINUE 
7  CONTINUE 
♦ 

OPENFILE  8,  "ALAT" 

REWIND  8 
♦ 

DO  9  J-1,4 

DO  10  1=1, 11 

READ  (8,  ) AALAT <  I,  J) 

10  CONTINUE 
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9  CONTINUE 
* 

OPENFILE  9,  HBLAT" 

REWIND  9 
* 

DO  11  J=l>  4 
DO  12  1  =  1,6 
READ  (9,  >BBLAT< I,  J) 

12  CONTINUE 
11  CONTINUE 
* 

CALL  MMULT  <  ALON,  AALON,  16,  6,  1,  V5) 
IF  (V.  LT.  5.  )  ALON <3,  l)=-32.  02 
IF  (V.  LT.  5.  )  ALON ( 7 ,  l)=-3.  37 

* 

CALL  MMULT (BLON,  BBLON,  9,  4,  1,  V3) 
BLON( 12,  1 ) =BLGN<  9,  1) 

BLON  <11, 1 )=BL0N<8,  1) 

BL0N(10,  1 ) =BLON<  7,  1) 

BLON  <7,  1 )  =0. 

BLON  <8,  1 )  =0. 

BLON  (9,  1  )=0. 

* 

CALL  MMULT  (ALAT,  AALAT,  11,  4,  1,  V3> 
DO  3  1=12, 16 
ALAT  <  1 ,  1  )=0. 

3  CONTINUE 
ALAT  <  14,  1)  =  1. 

* 

CALL  MMULT  (BLAT,  BBLAT,  6,  4,  1,  V3) 
BLAT (8,  1 ) =BLAT  <  6,  1) 

BLAT  <7,  1 )  =BLAT  ( 5,  1) 

BLAT (5, 1 ) =BLAT ( 4, 1) 

BLAT (4,  1 ) =BLAT ( 3,  1) 

BLAT  <3,  1  )=0. 

BLAT  <6,  1  )=0. 

DO  6  1=9, 12 
BLAT  ( 1 ,  1 )  =0. 

6  CONTINUE 
# 

RETURN 

END 


APPENDIX  F 

DETERMINATION  OF  STATE  FEEDBACK  GAINS 
("SFBCOF" ) 

This  subroutine  employs  a  table  look-up  algorithm  to 
determine  the  feedback  gains.  Feedback  gains  were  calculated 
for  three  design  points:  V  =  5,  30,  and  60  kt.  The  SFBCOF 
program  employs  the  following  logic: 


o< 

V< 

15 

USE 

V  =  5  DESIGN 

15* 

V< 

45 

USE 

V  =  30  DESIGN 

45  <  V 


USE 


V  =  60  DESIGN 


I 

I 

I 

I 

I 

I 

I 

I 

I 

J 

1 

I 

1 

1 

1 

I 

I 

i 

I 
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SUBROUT I ME  SFBCOF <  V , X KLON ,  X KL AT ) 

*  PROGRAMMER:  R.  V.  WALTERS 

*  PROGRAM  NAME:  SFBCOF 

*  SYSTEM:  DTSS  FORTRAN 

*  FUNCTION:  DETERMINES  AV-8B  HARRIER  SFB 

*  GAINS  AT  ANY  VELOCITY 

«• 

*********-*•«■■**•«■  ***■*•********■#■*■»••»■•*  *••«■**■#■**•*•** 

* 

* 

^LIBRARY  "KL0N5“ » "KL0N30", "KLONSO" 

•» 

DIMENSION  XKLON <  8 ) ,  XKLAT ( 8 ) 

* 

IF<  V.  GT.  15  )  GO  TO  2 

# 

ORENFILE  1,  " KL0N5 " 

REWIND  1 

OPENFILE  2,  "KLAT5" 

REWIND  2 
# 

DO  1  1  =  1,  8 
READ ( 1 , ) XKLON < I ) 

READ (2, ) XKLAT ( I ) 

1  CONTINUE 
RETURN 

■* 

2  CONTINUE 
•* 

IF  <  V.  GT.  45.  )  GO  TO  4 

* 

OPENFILE  3,  " KL0N30 " 

REWIND  3 

OPENFILE  4  "KLAT30" 

REWIND  4 
•* 

DO  3  1=1, 8 
READ <3, ) XKLON ( I ) 

READ (4, ) XKLAT < I ) 

3  CONTINUE 
RETURN 

« 

4  CONTINUE 

* 

OPENFILE  5, "KL0N60" 

REWIND  5 

OPENFILE  6, "KLAT60" 

REWIND  6 
* 

DO  5  1  =  1,8 
READ <5, ) XKLON < I ) 
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READ <6,  ) XKLAT  < I ) 
CONTINUE 

RETURN 

END 
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APPENDIX  G 

"DIGISIM"  DOCUMENTATION 


DIGISIM,  programmed  by  Associate  Professor  E. 


E.  Mitchell 


is  the  digital  simulation  utility  program  on  which  my  DTSS  model 
was  based.  This  appendix  describes  its  features  and  user-defined 
options. 


DIGISIM 


05/05/7? 


16:  35:  45 


PAGF 
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******  DIGISIM  -  DIGITAL  SIMULATOR  —  E.  E. 

* 

***  REMEMBER  —  RUN  FORTRAN 
« 


*  LINES  400-650: 

*  LINES  750-1000: 

*  LINES  1200-1500: 

*  LINES  1700-1850: 


INITIAL  TIME  DATA  BLOCK 
SYSTEM  DYNAMIC  EQUATIONS 
OUTPUT  STATEMENTS 
TERMINAL  CALCULATION  REGION 


MITCHELL  2/7 


* 

*  BASIC  DATA:  TMAX, DELT, NSAMPL 

*  BASIC  STATEMENTS:  X=ENTGRL ( XD,  X ( 0 ) > ,  CALL  OUTPUT < X,  " ) 

*  ADDITIONAL  FUNCTIONS:  STEP, FRSTOR, SCNDOR, PULSE, SQUARE, DELAY ,  DIGTAL 

*  PTRAIN 

*  FOR  MORE  DETAILS  CONTINUE  LISTING 

* 

****  INITIAL  TIME  REGION 


* 


*  THIS  PROGRAM  HAS  AN  INITIAL  TIME  DATA  REGION  EXTENDING  FROM  LINES  400-7 

*  THIS  SPACE  IS  USED  TO  DEFINE  CONSTANTS,  INITIAL  VALUES  AND  OTHER  NO1, 

*  CHANGING  PARAMETERS  REQUIRED  BY  YOUR  PROGRAM.  ALSO  INITIAL  CALCULATIONS 

*  MAY  BE  MADE  HERE  (IF  YOU  DEFINE  THE  REQUIRED  PARAMETERS). 

*  IN  ADDITION,  SPECIAL  CONSTANTS  FOR  BUILT  IN  DIGISIM  FUNCTIONS  ARE  BEFR 

*  HERE.  FOR  INSTANCE,  COMMONLY  DEFINED  TERMS  IN  THIS  V EA  ARE 


* 

* 

* 


* 

* 

* 

* 

* 

* 

* 


HE AD= " 


ONE  LINE  OF  80  OR  LESS  CHARACTER  BETWEEN  rHE  QUOTE 


DELT 

TMAX 

NSAMPL 

ZK 

ZA( 1-6 
TSAM 
TLAG 
TO 

TKPLOT 


THIS  IS  TYPED  OUT  AS  A  HEADING  FOR  THE  RUN 
INTEGRATION  STEP  SIZE  (DELT=0.  05  DEFAULT? 

TOTAL  SIMULATED  RUN  TIME  <  TMAX-3  DC:  JLT s 
OUTPUT  SPACING-  ANSWERS  PRINTED  EVERY  DEL  T  « > 4SAMPI 

TIME  INTERVAL  <NSAMPL=2  DEFAULT) 

DIGITAL  CONTROLLER  GAIN 

) »  ZB(  1-6)  DIGITAL  CONTROLLER  COEFFICIENTS  (  7.A<  ) --ZB<  >=0  DEFAULT 

DIGITAL  CONTROLLER  SAMP  TIME  (TSAM  =0  DEFAULT) 

DELAY  PERIOD  FOR  TIME  DELAY  ( TLAG=0  DEFAULT) 

INITIAL  TIME  IF  OTHER  THAN  0  (T0=0  DEFAULT) 

IF  TKPL0T=1  THE  OUTPUT  RESPONSE  TS  WRITTEN  INTO  A  TILE 


* 

* 

* 

* 

* 

**** 


SAVED  BY  YOU.  IT  IS  IN  THE  '-PER  FORM  F  OR  PLOTTING 
ON  THE  TEKTRONIC  GRAPHIC  TERMINALS.  ( l KPL0T=0  DEFAUL 
PPLANE  IF  PPLANE=1  THE  TKPLOT  FILE  IS  WRITTEN  TO  MAKE  A  PHASE 

PLANE  PLOT  OF  THE  FIRST  TWO  CALL  OUTPUT  VARIABLES 

DYNAMIC  REGION 


* 


*  THE  PROGRAM  HAS  A  DYNAMIC  REGION  FROM  STATEMENT  NUMBERS  750-1000.  IN  Tl 

*  AREA  THE  SYSTEM  DYNAMICS,  CONSISTING  OF  DIFFERENTIAL  ■  ^  iL  ’  BRAT  7 C 

*  EQUATIONS,  ARE  DEFINED. 

*  IN  THE  DYNAMIC  REGION,  THE  FULL  REGULAR  FORTRAN  LIBRARY  I S  'VARIABLE. 

*  FOR  THAT  MATTER,  IT  IS  AVAL  I ABLE  ANYWHERE  IN  THE  PROGRAM  IN  ADDITION 

*  THE  FOLLOWING  FUNCTIONS  ARE  FURNISHED:  ( X* INPUT, Y=0UTPU7 

* 

*  T  =  TIME 

*  Y  =  ENTGRL<  X, XO)  —  XO=INITIAL  VALUE  -  DECIMAL  POINT  REQUIRED  IF  NUN, 
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« 

« 

* 


« 

« 

* 

* 

* 

* 

« 

* 

* 

* 

* 

* 

* 

* 

* 


Y  =  STEP  <  T1 )  —  Y=0  FOR  T<T1,  Y=1  FOR  T>=T1 

Y  =  PULSE <  T 1 »  T2 )  —  Y=0  EXCEPT  Y=1  WHEN  Tl<=TOT2 

Y  =  SQUARE ( T 1 »  T2»  T3»  P )  —  SQUARE  WAVE  WITH  AMPLITUDE  OF  Y=l.  LEADING 

OF  FIRST  PULSE  AT  T3,  PULSE  WIDTH  IS  T1  TO 
REPEATS  EVERY  P  SECONDS. 

Y  =  FRSTOR(X, A, B, C, D)  --  FIRST  ORDER  TRANSFER  FUNCTION  DEFINED  AS 

<A*S+B)/<C*S+D)  C  IS  NOT  ZERO 

Y  =  SCNDOR<X, A, B, C, D, E, F)  —  SECOND  ORDER  TRANSFER  FUNCTION  DEFINED  ( 

<A*S*S+B*S+C)/(D*S*S+E*S+F)  D  IS  NOT  ZEf 

Y  =  DELAY ( X,  TLAG )  —  TIME  DELAY,  Y  =  X  AFTER  X  IS  DELAYED  TLAG  SECONI 

Y  =  DIGTAL(X)  —  DIGTAL  FILTER  OR  CONTROLLER  SIMULATOR.  DELT  SHOULD  L 

OTSAM/IO.  TSAM=SAMPLE  INPUT-OUTPUT  TIMES  OF  THE  DI( 
CONTROLLER.  TRANS  FUNCTION  HAS  6  ZA,  ZB  TERMS  AND  IS 
<  ZK*<  ZA  < 1 ) *X+ZA( 2 ) *X-1+ZA<  3) *X-3)  -  ZB<2)*Y-1 
-ZB( 3) *Y-2  )/ZB(l) 

X=CURRENT  INPUT,  X-l  =  LAST  INPUT,  ETC.  DITTO  FOR  Y 

Y  =  LIMIT(U1, X, U2)  —  Y=U1  IF  XCUl,  Y=U2  IF  X>U2,  Y=X  OTHERWISE 

Y  =  PTRAIN(Tl)  --  Y=1  IF  T  =  N*T 1 »  THIS  IS  A  PULSE  TRAIN,  UNIT  PULSE 

EVERY  N*T1  SECONDS 


* 

* 


OUTPUT  REGION 


* 


*  THE  PROGRAM  HAS  AN  OUTPUT  REGION  FROM  STATEMENT  NUMBER  1200-1500.  IN  Th 

*  REGION  THE  VARIABLES  TO  BE  PRINTED  AND  PLOTTED  ARE  DEFINED.  TO  OUTPUT 

*  X,  XDOT  AND  X2D0T ,  THE  STATEMENTS  ARE. 

* 

*  CALL  OUTPUT (X,  "X") 

*  CALL  OUTPUT<  XDOT,  "XDOT") 

*  CALL  OUTPUT (X2D0T, "ACCEL" ) 


*  ANY  8  CHARACTERS  MAY  BE  USED  BETWEEN  THE  "  ",  THESE  ARE  THE  COLUMN 

*  HEADINGS  OVER  THE  PRINTED  VALUES. 

**  NOTE  ** * ****  TIME  IS  OUTPUT  AUTOMATICALLY 

* 

**  NOTE  **  THE  FIRST  6  CALL  OUTPUT  VARIABLES  ARE  PLOTTED 
* 

****  TERMINAL  CALCULATION  REGION 
* 

*  THE  PROGRAM  HAS  A  TERMINAL  REGION  FROM  STATEMENT  NUMBERS  1700-1850.  TH 

*  AREA  IS  USED  TO  MODIFY  PARAMETERS.  ECT.  AT  THE  END  OF  A  RUN.  ONE  CAN 

*  THEN  GO  BACK  AND  INTEGRATE  AGAIN  AND  AGAIN.  A  TYPICAL  PROBLEM  THAT 

*  WOULD  USE  THIS  AREA  IS  A  BOUNDARY  VALUE  PROBLEM. 

* 

*  IF  YOUR  PROGRAM  CONTAINS  THE  STATEMENT 

*  GO  TO  9000 

*  CONTROL  IS  TRANSFERED  BACK  TO  THE  TOP  OF  THE  INITIAL  AREA 

*  IN  A  LIKE  MANNER,  THE  STATEMENT 

*  GO  TO  9100 

*  TAKES  CONTROL  TO  THE  BOTTOM  OF  THE  INITIAL  AREA, (TOP  OF  DYNAMIC  AREA) 

* 
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APPENDIX  H 

"AHEAD"  -  ADVANCED  HARRIER 
ELECTRONIC  AUGMENTATION  DEVICE 

AHEAD  is  the  DIGISIM-based  digital  Harrier  simulation.  he 
output  traces  of  Appendix  L  were  made  using  AHEAD  in  conjunction 
with  a  special  APL  graphing  routine.  AHEAD  is  an  extremely  versatile 
design  tool.  The  user  must  program  control  inputs  and  initial 
conditions  prior  to  run  time.  At  run  time,  user  inputs  TMAX,  his 
desired  SAS  option,  desired  output  format,  and  output  file.  System 
response  data  is  then  output  to  the  user-defined  output  file. 
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*  PROGRAMMER:  R.  V.  WALTERS 

*  PROGRAM  NAME:  AHEAD 

*  SYSTEM:  DTSS  FORTRAN 

*  FUNCTION:  SIMULATES  THE  LAT  AND  LONG 

*  DYNAMICS  OF  THE  AV-8B  HARRIER, 

*  INCLUDING  BOTH  EXISTING  STABILITY 

*  AUG  AND  A  STATE  FB  CONTROLLER 

* 

* 

w****************************************** 

LIBRARY  "MATCOF " 

LIBRARY  "SFBCOF " 

IMPLICIT  REAL< I-N) 

INTEGER  ICZ, NICZ,  ITZ, NITZ,  ITMAX,  KERR,  IAZ, NZZ, KPZ, NSAMPL, KTZ,  IIZ,  JJZ,  IZ 

INTEGER  ICNT , NVAR,  IDIGL,  ICN,  IADCNT,  ISZZ 

DIMENSION  ALON<  16,  1 ) ,  ALAT  (16,  1),BL0N<12,  1),BLAT<12,  1) 

COMMON  ZGUT, XVZ<6, 402),  ICNT,  NVAR, T,  TO,  TMAX,  DELT,  HDELT ,  ICZ,  NICZ,  ITZ,  NITZ, 

&  ITMAX, KERR,  IAZ( 100),  Y1PZ ( 100) ,  Y4PZ < 100 ) ,  X1PZ( 100),  X4PZ < 100) ,  NZZ,  KPZ, 

&  NSAMPL,  KTZ 

DIMENSION  XKL0N<8), XKLAT<8) 

COMMON  /HOLD/ 1 ADCNT ,  ITQA, JCNT,  I QCNT , HLD  <  500 ) 

COMMON  /0UTP2/  IIZ, JJZ, AZZ,  ISZZ,  AMXZ < 25) ,  AMNZ ( 25) ,  HEAD,  NPLOT ,  TKPLOT ,  PPLANI 
COMMON/DIGITA/ IDIGL,  ZA<6),  ZB<  6) ,  TSAM,  ZK,  NCNT 

* 

CHARACTER  YES 

CHARACTER  PDQ, PDQ2,  IPQ.  AZZ < 25 ) ,  HEAD*80 

REAL  LIMIT 

INTEGER  TKPLOT 

DATA  ICN,  IIZ,  JJZ/O,  0,  0/ 

Z0UT=ICNT=IDIGL=0 
TKPL0T=0 
NPL0T=6 
DO  109  IZ=1, 4 
ZA  (  I Z  )  =0. 

109  ZB(IZ)=0. 

ZK=0. 

I ADCNT=0 
DO  108  IZ=1, 500 
108  HLD(IZ)=0. 

1  DO  2  IZ-1, 100 

2  I AZ  < I Z ) =0 
9000  CONTINUE 

to»0. 

DELT=.  05 

NICZ«10 

NITZ=4 

*###  INITIAL  REGION  *** 

* 

PRINT,  "IC'S  *<  INPUTS  CORRECT?" 

PRINT, 

PRINT,  "IF  NOT, STOP!" 
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PRINT, 

PRINT, 

* 

PRINT, "INPUT  TMAX " 

INPUT, TMAX 
PRINT, 

DELT*=.  002 
NSAMPL=50 

HEAD= "  AV-8B  TIME  RESPONSE" 

TKPL0T=1 

VMUL= 1 . 

* 

XL0N10=0. 

XL0N20=0. 

XL0N30=0. 

XL0N40=0. 

* 

XLAT10-0. 

XLAT20-0. 

XLAT30=0. 

XLAT40=0. 

* 

QSAS=0. 

RSAS=0. 

PSAS=0. 

SFBL0N=0 

SFBLAT-O. 

* 

PRINT, "CONVENTIONAL  SAS  ON?" 

* 

INPUT,  YES 

IF  < YES.  NE.  "YES")  GO  TO  506 
QSAS= 1 . 

RSAS=1. 

PSAS=1. 

* 

506  PRINT, " SFB  SAS  ON?" 

INPUT,  YES 

IF  (YES.  NE.  "YES")  GO  TO  516 
SFBL0N=1. 

SFBLAT-1. 

* 

516  CONTINUE 

* 

IFLAG1=1 

IFLAG2-1 

PRINT, "SUPPRESS  LONG  DYNAMICS?" 
INPUT, YES 

IF  (YES.  EQ.  "YES"  )  I  FLAG  1=0 
PRINT,  "SUPPRESS  LAT  DYNAMICS?" 
INPUT,  YES 
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IF  (YES.  EQ.  "YES"  )  IFLAG2=0 
PRINT, 

* 

PRINT, "INPUT  VELOCITY" 

INPUT,  V 
VINIT=V 
PRINT, 

* 

CALL  MATCQF  <  V, ALON,  BLON,  ALAT,  BLAT) 

CALL  SFBCOF(V, XKLON,  XKLAT) 

* 

PR I NT, " I N I T I AL  COEFF I C I ENTS  CALCULATED " 
PRINT, 

* 

PRINT, "BEGIN  SIMULATION" 

PRINT, 

•* 

*  LET  T=0? 

* 

* 

* 

* 

*«#« 

IF ( TKPLOT  .  EQ.  0)G0  TO  9100 
WRITE <  ZOUT , S7 ) 

87  FORMAT <"  WHAT  IS  YOUR  SAVED  FILE  NAME") 

READ <  ZOUT,  14 ) IPQ 

OPENFILE  3,  IPQ 

REWIND  3 

ENDFILE  3 

9100  KTZ=0 

T=TO 

HDELT=0.  5*DELT 

KERR=0 

ICZ=1 

ITMAX=NICZ 

4  I TZ= 1 

5  NZZ=0 

****  DYNAMIC  REGION  *** 

IF(  IFLAG1.  EQ.  0)  GO  TO  838 
* 

*************************** 

*  LONGITUDINAL  SIMULATION 
*************************** 
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* 

BKL0N1=XKL0N< 1 ) *XL0N1 +XKL0N<  2 ) *XL0N2+XKLQN  <  3 ) *XL0N3+XKL0N<  4 ) *XL0N4 
BKL0N2=XKL0N<  5 )*XL0N1 +XKL0N<  6 ) *XL0N2+XKL0N<  7 ) *XL0N3+XKL0N ( 8 ) *XL0N4 
* 

762  THR0T*0  ! THROTTLE 
UL0N3®0.  • JET  ANGLE 

* 

UL0N2»THR0T-SFBL0N*BKL0N2 

* 

LONSTK-O. 

L0NSTK»LIMIT<-8.  ,  LONSTK,  4.  )  !  LONSTK 

* 

Q1DD-17410.  *XL0N4-132.  *Q1D-17410.  *Q1 
Q1D»ENTGRL<Q1DD,  O.  ) 

Q1=ENTGRL<Q1D,  0.  ) 

# 

BXLOND-4.  04*(QlD+5.  3*Q1)-1.  35*BXL0N 
BXLON=ENTGRL ( BXLOND, O.  ) 

BXLON»QSAS*LIMIT < 84,  BXLON,  .  84) 

* 

ELON*LONSTK-BX  LON-SFBLON*BKLON 1 
* 

ULDN1D-12.  *  <  EL0N-UL0N1  ) 

UL0N1*ENTGRL(UL0N1D,  0.  )  !  ELEV/STK  INPUT 

UL0Nl»LIMIT<-8.  84,  UL0N1,  4.  84) 

* 

XL0N1D=AL0N<1,  1 ) *XL0N1+ALQN< 2,  1 ) *XL0N2+AL0N< 3,  1 )*XLQN3+AL0N< 4,  1)*XL0N4 
XL0N1D=XL0N1D+BL0N< 1,  1 )*UL0N1+BL0N<2,  1 > *UL0N2+BL0N< 3,  1)*UL0N3 
XL0N2D*AL0N  <  5,  1 ) *XL0N1+AL0N< 6,  1 > *XL0N2+AL0N< 7,  1 >*XL0N3+AL0N<8,  1)*XL0N4 
XL0N2D*XL0N2D+BL0N ( 4,  1 >*UL0N1+BL0N<5,  1 ) *UL0N2+BL0N< 6,  1)*UL0N3 
X L0N3D=AL0N < 9,  1 ) *XL0N1+AL0N< 10,  1 >*XL0N2+AL0N< 1 1 ,  1 ) *XL0N3+AL0N<  12,  l)*XL0i>14 
XL0N3D=XL0N3D+BL0N<7,  1 )*UL0N1+BL0N<8,  1 )*UL0N2+BL0N<9,  1)*ULGN3 
XL0N4D=AL0N< 13,  1 )*XL0N1+AL0N< 14,  1 )*XL0N2+AL0N< 15,  1 )*XL0N3+AL0N< 16,  1)*XLCN' 
XL0N4D=XL0N4D+BL0N< 10,  1 ) *UL0N1+BL0N< 11,  1 ) *UL0N2+BL0N< 12,  1 )*UL0N3 
* 

XL0N1«ENTGRL( XL0N1D,  XL0N10) 

XL0N2=ENTGRL ( XL0N2D,  XL0N20) 

XL0N3*ENTGRL < XL0N3D,  XL0N30) 

XL0N4»ENTGRL<XL0N4D,  XL0N40) 

* 

838  CONTINUE 

IF ( IFLAG2.  EQ.  0)  GO  TO  978 
* 

*«•»*****#*•»**•**#***** 

*  LATERAL  SIMULATION 


BKLAT1-XKLAT ( t ) *XLAT1+XKLAT <  2 ) *XLAT2+XKLAT  <  3 ) *XLAT3+XKLAT  <  4 ) *XLAT4 
BKLAT2-XKLAT  <  5 ) *XLAT1 +XKLAT  <  6 ) *XLAT2+XKLAT  <  7 ) *XLAT3+XKLAT  <  8 ) *XLAT4 
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ULAT3=0.  !  NO  INPUT 

* 

P1DD=17410.  *XLAT2-132.  *P1D-17410.  *P1 
P1D=ENTGRL<P1DD,  0.  ) 

P1=ENTGRL<P1D,  0.  ) 

* 

P2D=4.  07*<P1D+1.  28*P1)-.  774*P2 
P2=ENTGRL<P2D,  0.  ) 

P2=PSAS*LIMIT<-  68 j P2,  .  68) 

* 

LATSTK=0.  !  LATSTK 

LATSTK=LIMIT <  -4.  ,  LATSTK,  4.  ) 

* 

EPLAT=LATSTK-P2-SFBLAT*BKLAT 1 
* 

ULAT1D=25.  * ( EPLAT-ULAT 1 ) 

ULAT1=ENTGRL ( ULATIDi  0.  )  !  ALRN/STK  INPUT 

UL AT  1=LIMIT(— 4.  68,  ULAT1,4.  68) 

•* 

R1DD=17410.  *XLAT3— 132.  *R  ID- 174 10.  *R1 
R1D=ENTGRL(R1DD,  0.  ) 

R1=ENT6RL(R1D,  0.  ) 

* 

R2D=23.  *  ( R 1  — R2 ) 

R2=ENTGRL  <  R2D,  0.  ) 

* 

R3D=1  73*R2D-  32*R3 
R3=ENTGRL  ( R3D,  0.  ) 

* 

R4D=23.  26*  <.  24*ULAT1-R4) 

R4=ENTGRL  <  R4D,  0.  ) 

* 

AY=— 1.  *XLAT  1D+2.  32*XLAT2D-1 1.  *XLAT3D+32.  2*XLAT4 
* 

AY1DD=318.  *AY— 19.  45*AY1D-164.  *AY1 
AY  1  D=ENTGRL  ( AY  1 DD,  0.  ) 

AY  1  =ENT6RL ( AY  1 D,  0.  ) 

* 

AY2D=.  12*<  AY1D+4.  *AYl)-8.  *AY2 
AY2=ENTGRL<  AY2D,  0.  ) 

* 

BRLAT =R4-R3-AY2 

BRLAT =RSAS*L I M I T < -.  7, BRLAT,  .  7) 

* 

RUDPED=0.  RUDPED 

RUDPED=LIMIT  <  -2.  ,  RUDPED,  2.  ) 

* 

ERLAT-RUDPED+BRLAT-SFBLAT*BKLAT2 

* 

ULAT2D=57.  *  <  ERLAT-ULAT2 ) 

ULAT2=ENTGRL < ULAT2D, 0.  )  !  RUDR/STK  INPUT 
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ULAT2»LIMIT<-2.  7,  ULAT2,  2.  7) 

* 

XLAT1D«ALAT<  1,  1  )*XLATH*LAT<2,  1 )*XLAT2«*LAT<3,  1 )  *XLAT3+ALAT  <  4,  1)*XLAT4 
XLAT1D*XLAT1D+BLAT<1,  1 )*ULAT1+BLAT<2, 1 ) *ULAT2+BLAT  <  3, 1)*ULAT3 
XLAT2D-ALAT<5,  1  )*XLATDALAT<6,  1  )*XLAT2«-ALAT<7,  1  )*XLAT3+ALAT<8,  1>*XLAT4 
XLAT2D-XLAT2IH-BLAT  <  4,  1 )  *ULAT  1 +BLAT  <  5#  1 )  *ULAT2+BLAT  <  6,  1 )  *ULAT3 
XLAT3D*ALAT<9.  1 ) *XLAT1+ALAT< 10, 1 ) *XLAT2+ALAT  <11.  1 )*XLAT3+ALAT< 12.  1)*XLAT4 
XLAT3D«XLAT3D«-BLAT<7,  1  >*ULAT1+BLAT<8,  1 >*ULAT2+BLAT<9,  1)*ULAT3 
XLAT4D-ALAT!  13.  1  >*XLAT1+ALAT<  14,  1 )*XLAT2+ALAT<  15,  1  )*XLAT3+ALAT  <16,  1  >#XLAT* 
XLAT4D«XLAT4D-HBLAT<  10,  1 )  *ULAT  1 +BLAT  <11,  1 )  *ULAT2+BLAT  ( 12,  1  )*ULAT3 
* 

XLAT1-ENTGRLCXLAT1D,  XLAT10) 

XLAT2-ENTGRL  <  XLAT2D,  XLAT20 ) 

XLAT3*ENTQRL  <  XLAT3D,  XLAT30 ) 

XLAT4-ENTGRL ( XLAT4D,  XLAT40 ) 

• 

978  CONTINUE 
* 

IF  <  T.  NE.  1.  )  00  TO  1003 

******** 

wwwW 

IF<NF1.  NE.  DPRINT,  "T«l" 

NF1»1 

1003  IF  IT.  NE.  2.  )  00  TO  1006 
IF<NF2.  NE.  DPRINT,  "T«2" 

NF2*1 

1006  IF  <T.  NE.  3.  )  00  TO  1009 
IF (NF3.  NE.  DPRINT,  "T«3" 

NF3»1 

1009  IF  <  T.  NE.  4.  )  GO  TO  1012 
IF ( NF4.  NE.  1  > PRINT,  "T-4" 

NF4*1 

1012  IF  <T.  NE.  TMAX)  GO  TO  1015 
IF (NF5.  NE.  DPRINT,  "T=THAX" 

NF5-1 

1015  CONTINUE 

IF  (KERR  .  NE.  0)  00  TO  9 

ITZ-ITZ+1 

IF < ITZ  .  L€.  I TMAX )  00  TO  5 
KPZ-0 

***  OUTPUT  REGION  **#  LINES  1200-1500 

« 

IFdFLAOl.  NE.  1)00  TO  1245 
* 

CALL  OUTPUT <  XL0N1,  "U<  X-VEL) " ) 

CALL  OUTPUT <  XL0N2,  "W<  Z-VEL ) " ) 

CALL  OUTPUT< XL0N3,  "THETA") 

* 

IF <  IFLAG2.  EQ.  1)00  TO  1245 
* 

CALL  OUTPUT! XL ON 4,  "Q<DPIT)M> 
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1245  CONTINUE 

*  *  * 

IF <  IFLAG2.  NE.  1)  GO  TO  1267 

CALL  OUTPUT  <  XLAT 1 ,  " V  < Y-VEL ) " ) 

* 

IF<  IFLAG1.  EQ.  1)  GO  TO  1263 
* 

CALL  OUTPUT  <  XLAT2*  "P  < DROLL ) " ) 

« 

1263  CONTINUE 

CALL  OUTPUT ( XLAT3/  " R ( DYAW ) " ) 

CALL  OUTPUT <XLAT4,  "PHI" ) 

* 

1267  CONTINUE 
* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

**** 

JJZ=JJZ+2 

IF ( ICZ  .  EQ.  0)GG  TO  6 
ICZ=0 

ITMAX=NITZ 

6  IF ( T  .  LT.  (TMAX+.  00001) > GO  TO  68 
IF ( TKPLOT  .  NE.  0)G0  TO  7 

WRITE ( ZOUT ,  12) 

12  FORMAT ("  WANT  TO  INCREASE  TMAX " ) 

READ<  ZOUT i  14) IPQ 

14  FORMAT  <  V ) 

IF < IPQ  .  NE  "YES")  GO  TO  7 
WRITE<  ZOUT/  13) 

13  FORM  AT  ("  TMAX35"  ) 

READ<  ZOUTi  14) TMAX 
GO  TO  6 

68  KTZ=KTZ+1 
T=KTZ 

T=T*DELT+TO 
GO  TO  4 

9  WRITE<  ZOUT »  10) 

10  FORMAT < 5X/  "  INTEGRATION  ASSIGNMENT  ERROR") 
STOP 

7  ISZZ-2 

»***  TERMINAL  REGION  **# 


I 
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* 

* 

* 

* 

#*#* 

I F  < TKPLOT  .  EQ.  0)CALL  RANGE 
IF  <  TKPLOT  .  EQ.  0)  CALL  PL0T2 
IF (TKPLOT  .  NE.  0)  CALL  TEKPLT 
STOP 
END 

FUNCT I ON  ENTGRL ( X ,  YO ) 

COMMON  OUT,  XVZ  ( 6,  402),  ICNT,  NVAR,  T,  TO,  TMAX,  DELT,  HDELT ,  IC,  NIC,  IT,  NIT,  ITMAX, 

MAI 100), Yl( 100), Y2( 100), XI ( 100), X2< 100), N, KP, NSAMPL, KT 

N=N+1 


IF  < 

IC. 

EQ.  0) 

GO 

TO 

20 

IF  < 

:  it. 

NE.  1) 

GO 

TO 

10 

IF  ( 

:  IA(N).  NE. 

0) 

GO 

TO 

50 

I A  ( 

:n)  = 

=  1 

GO 

TO 

15 

10 

IF  ( 

I A<  N) . 

EQ. 

0) 

GO 

TO 

50 

15 

Y2N=Y0 

GO 

TO 

40 

20 

IF  ( 

IAIN). 

EQ. 

0) 

GO 

TO 

50 

IF  * 

:  it. 

NE.  1) 

GO 

TO 

30 

*  PREDICTOR 
Y1N=Y2IN) 

X 1N=X2I N) 

Y1 (N)=Y1N 
XI <N)=X1N 
Y2N=Y 1 N+X 1 N*DELT 

*  CORRECTOR 
30  Y2N=Y1IN)+IX1IN)+X)*HDELT 
40  Y2IN)=Y2N 
X2<N)=X 
ENTGRL=Y2N 
RETURN 
50  KERR=1 
RETURN 
END 

FUNCTION  STEP ( T 1 ) 

COMMON  OUT, XVZ <6,  402),  ICNT, NVAR, T 
Y-l. 

IFIT.  LT.  Tl)  Y=0. 

STEP=Y 
RETURN 
END 

FUNCTION  PULSE <T1,T2) 

COMMON  OUT,  XVZ (6,  402) ,  ICNT,  NVAR,  T 
Y=0. 

IF < T.  GE.  Tl.  AND.  T.  LE.  T2)  Y-l. 
PUL.SE=Y 
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RETURN 

END 

FUNCTION  SQUARE < T 1, T2, T3) 

COMMON  OUT, XVZ  <  6, 402),  ICNT,  NVAR,  T,  TO,  TMAX,  DELT,  HDELT 
Y«0. 

TM=AMOD  <  T-T 1 ,  T3) 

IF(  TM  .  GE.  0.  .AND.  TM  .  LE.  T2)  Y  =  1. 

SQUARE=Y 

RETURN 

END 

FUNCTION  FRSTOR<  X,  A,  B,  C,  D) 

COMMON  /FRST/  SZ,  Z 
IF< C.  EQ.  ZERO)  GO  TO  1 
SZ=X-D*Z/C 
Z=ENTGRL  ( SZ ,  0.  ) 

FRSTOR= <  A*SZ+B*Z ) /C 
RETURN 

1  PRINT  2 

2  F0RMAT<5X,  "FIRST  ORDER  TRANSFER  FUNCTION  OUT  OF  ORDER") 

STOP 

END 

FUNCTION  SCNDOR(X,  A,  B,  C,  D,  E,  F) 

COMMON  /SCNB/  S2Z, SZ,  Z 
IF< D.  EQ.  ZERO)  GO  TO  1 
S2Z=X-(E*SZ+F*Z)/D 
SZ=ENTGRL<S2Z,  0.  ) 

Z=ENTGRL<SZ,  0.  ) 

SCNDOR= ( A*S2Z+B*SZ+C*Z ) /D 
RETURN 
1  PRINT  2 

?  FORMAT <5X,  "SECOND  ORDER  TRANSFER  FUNCTION  OUT  OF  ORDER") 

STOP 

END 

SUBROUT I NE  OUTPUT  <  X , A 1 ) 

CHARACTER  A 1 , A  <  25 ) ,  HE AD*80 

COMMON  OUT,  XVZ  (6,  402),  ICNT,  NVAR,  T,  TO,  TMAX,  DELT,  HDELT,  IC,  NIC,  IT,  NIT,  I  TMAX, 
&I A( 100), Y1 < 100),  Y2< 100),  XI < 100),  X2< 100),  N,  KP,  NSAMPL,  KT 
COMMON  /0UTP2/I I , JJ,  A,  ISZZ, AMXZ < 25 ) ,  AMNZ < 25 ) ,  HEAD,  NPLOT ,  TKPLOT,  PPLANE 
COMMON  /0UTP1/  V(25),  IP, KASE, KS, TSAV,  ICN 
* 

IF< T  .  NE.  TO) GO  TO  10 
KS  =  0 

IF ( JJ  .  GT.  0)  GO  TO  40 
11=11+1 
A( I I )=A1 
ISZZ=1 

NVAR  =  MIN ( II,  6,  NPLOT) 

GO  TO  40 

10  IF< KP  .  NE.  0)  GO  TO  40 
IF<KS  .  EQ.  0)  GO  TO  60 
20  KS=KS+1 


" 
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IF (KS  .  GE.  NSAMPL)  KS  =  0 
40  KP=KP+1 

I F  <  K'S  .  NE.  0)  RETURN 
TSAV  ■  T 
V<KP)  =  X 
RETURN 

*  PRINT  LAST  TIME  PERIOD 
60  CONTINUE 
1CNT=ICNT+1 

DO  120  1=1,  N V AR 

120  XVZ(I, ICNT )  =  V(I> 

IF ( TKPLOT  .  NE.  0)G0  TO  20 
IF <  ICN  .  EQ.  0)  WR I TE < OUT ,  900 ) 

I F <  I CN  .  EQ.  0)  PRINT,  HEAD 
IF <  ICN  .  LT.  51 )  GO  TO  70 
ICN  =  0 

WRITE<OUT> 915) 

WRITE (OUT,  900) 

70  CONTINUE 
KT1=MIN(5, II) 

I F  ( J  J  EQ.  2)  WRITE  (OUT,  910)(A<I),  1  =  1,  KT1 ) 
IF( 1 1  .  GT.  5)WRITE(0UT,  905) 

WRITE (OUT, 920) ICNT, TSAV, (V(I), 1=1, KT1) 

I CN= I CN+ 1 

I F  ( 1 1  .  LE.  5)  GO  TO  100 

* 

ICN  =  ICN  +  1 
80KT1=KT 1  +  1 
KT2=MIN(KTl+4,  II) 

I F  ( J  J  .  EQ.  2 )  WR  I  TE  ( OUT ,  930 )  ( A  ( I ) ,  I  =KT  1 ,  KT2 ) 
WR I TE ( OUT , 940) ( V( I ) ,  I=KT1, KT2) 

I CN= I CN+ 1 
KT1  =  KT2 

IF(  II  .  GT.  KT2)  GO  TO  80 
100  JJ=4 
CALL  RANGE 
GO  TO  20 

900  FORMAT ( 1 H 1 ,////// > 

905  FORMAT (1H  ) 

910  FORMAT (/I IX,  "TIME",  3X,  5(2X,  A8,  2X ) ) 

915  FORMAT <////> 

920  FORMAT ( 1H  ,  13,  IX,  6G12.  5) 

930  FORMAT ( /,  18X, 5(2X,  A8,  2X ) ) 

940  FORMAT ( 1H  ,  15X,  5G12.  5) 

END 

SUBROUTINE  PL0T2 
* 

* 

*  N=  NUMBER  OF  POINTS 

*  NVAR=  NUMBER  OF  PLOTS 


10 
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* 

* 

*  V=  NVAR  X  N  MATRIX;  EACH  COLUMN  CORRESPONDS  TO  A  VECTOR 

*  COMPOSED  OF  POINTS  TO  BE  PLOTTED 

* 

* 

*  MAX  NVAR  =  6 

*  MAX  N  =  400  WITHOUT  CHANGING  DIMENSIONS 

* 

* 


COMMON  OUT, V ( 6, 402 ) , N, NVAR 

COMMON  /0UTP2/I I ,  JJ, A2,  ISZZ,  AMXZ < 25 ) ,  AMNZ ( 25) ,  HEAD,  NPLOT ,  TKPLOT. PPLANE 
DIMENSION  SF (7), IS(7),H(7), IC(7) 

CHARACTER  BLANK,  DOT,  STR,  MINS, HEAD*80 
CHARACTER  X ( 8) ,  LINE ( 120) ,  A2< 25 ) 

BLANK="  " 

DOT=" I " 

STR="+" 

MINS="-" 

X<  1  )=,•  +  ,, 

X  ( 2 )  =  "  *  " 

DO  9999  1  =  1,  7,  2 
9999  X(I)=X(1) 

DO  9998  1=2, 8, 2 
9998  X  < I ) =X  <  2 ) 

• 

* 

* 

DO  1  1=1, NVAR 
H< I ) =0 
SF  ( I )  =0.  0 
IS<I)=I 
1  IC(I)=1 

* 

DO  6  1=1, NVAR 

V(I,  1)=V(I,  1)  +  1.  E-6 
DO  6  J= 1 , N 
IF(VdfJ)  )3»  3,  4 

3  IC( I ) =2 

4  I F  <  ABS  ( V  ( I ,  J )  )  -H  ( I )  )  6,  6,  5 

5  H ( I ) =ABS  <  V ( I , J ) ) 

6  CONTINUE 

* 

* 

*  GENERATE  SCALE  FACTORS 

* 

* 

DO  7  1=1, NVAR 
A=&0/ ( NVAR* I C ( I ) ) 

SF  <  I )  =A/  <  H  ( I )  +.  1*H(I)  ) 


AHEAr 


A=(30  *( IC( I )-l )+60  * ( I— 1 ) ) /NVAR 
DO  7  J=l, N 

V(I, J)=SF(I)*V(I, J)+A 


WRITE  HEADING 


WR I TE ( OUT ,  500 ) NVAR 
IF ( NVAR— 4 ) 8,  8,  9 

8  NS=  1 
NQ=  1 
NP=NVAR 
GO  TO  10 

9  NS=2 
NQ=  1 
NP=4 

10  DO  11  K=l,  NS 

WR I TE ( OUT » 501 ) ( A2( I ) 
WR I TE  <  OUT , 502) (SF( I ), 
NQ=5 

11  NP=NVAR 
WRITE (OUT,  505) 

PRINT,  HEAD 

WRITE (OUT, 503) 


MAKE  AXIS 


LLNTH  =  61 

K=0 

DO  12  1=1, LLNTH 
LINE( I )=MINS 

WRITE (OUT, 504 ) K,  <LINE< I ), 1=1, LLNTH) 


BLANK  THE  LINE 


DO  13  1=1, LLNTH 
13  LINE( I ) =BLANK 


PLOT  THE  VARIABLES 


DO  15  .J=l,  N 
LINE( 1 )=DOT 
DO  14  1  =  1,  NVAR 
JK= ( 60  * I ) /NVAR  +1 
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JP=JK- <30  *  < I C  < I ) - 1 ) ) /NVAR 
JL=V  <  I ,  -J )  +2.  5 
L I NE ( JK ) =D0T 
LINE< JP)=D0T 

LINE( JL)=X< I ) 

WR I TE  <  OUT /  504 ) J,  <  L I NE ( K ) ,  K- 1 , LLNTH ) 
DO  15  K=l, LLNTH 
L I NE ( K ) =BL ANK 


RETURN 


500  FORMAT (  ////5X,  "SYSTEM  PLOT  IN"/5X,  12,  5X  "UfiRTm 

501  F0RMAT(5X.  -VARIABLES.  ",  5X,  4< A8,  6X»  VARIABLES  V) 

502  FORMAT (3X, "SCALE  FACTORS", 4<G14.  5) ) 

503  FORMAT <///lX,  6<  10X,  "1")) 

504  FORMAT (  IX.  14,  6X,  121A1 ) 

505  FORMAT <//) 

END 

FUNCTION  DELAY < C,  TLAG) 

CHARACTER  A  <  25 ) ,  HEAD*80 

COMMON  OUT,  XVZ  <6,  402) ,  ICNT,  NVAR,  T,  TO,  TMAX,  DELT, HDELT  If  WTr 

imanpl',  kt”AX'  KERR' Ifl(  100)'  V1  ‘ 100’'  Y4I100)' *>>  >t>o, ,  x4uoo>,  n.  kp, 

cotX!  "i it  HEfi0,  NRL0T' TKRL0T’ RPLflNE 

SS7,  "“Hi1 JCNT- ,QCNT- HLDI500) 

IADCNT=1 

ITQA=1 

•JCNT=0 

IQCNT=0 

102  IF < T.  GT.  TLAG)  GO  TO  103 
DELAY=HLD ( ITQA) 

HLD< ITQA)=C 
I QCNT = I QCNT + 1 
GO  TO  104 

103  I F  <  JCNT.  GE.  1)  GO  TO  106 
I TQA=1 

JCNT=1 
GO  TO  107 

106  IF  ( ITQA.  LT.  I  QCNT )  GO  TO  107 
JCNT=0 

107  DELAY=HLD< ITQA) 

HLD< ITQA)=C 

104  ITQA=ITQA+1 
RETURN 

END 


-133- 


AHEAD 


05/05/79 


20. 22. 18 


PAGE 


14 


FUNCTION  DIGTAL(Ul) 

COMMON  OUT, XVZ  <  6, 402),  ICNT, NVAR, T, TO, TMAX, DELT, HDELT ,  IC, NIC,  IT, NIT,  ITMAX, 
&IA< 100), Y1PZ  < 100), Y4PZ  < 100),  X1PZ( 100),  X4PZ< 100), N, KP, NSAMPL, KT 
COMMON/DI G I T A/ 1 D I GL , ZA<6), ZB<6>, TSAM,  ZK, NCNT 
COMMON/  EXTRA/  NTERM, U < 6 ) »  Y <  6 ) 

* 

if 

IF ( IDIGL  .  GT.  0)  GO  TO  100 
I D I GL= 1 
NCNT =0 
NTfcRM=6 

DO  10  1  =  1,  NTERM 
U(I)  =  O. 

10  Y(I)  =  0. 

«• 

•* 

*  ADJUST  DELT  TO  BE  INTEGRAL  MULTIPLE  OF  TSAM 
•* 

IF  (DELT  .  GT.  TSAM/10.  )  GO  TO  50 
PP=1. 

20  RR=TSAM/PP/DELT 
JJ=RR 

I P  (  <  RR-FLOAT  ( J J )  )  .  LT.  .001)  GO  TO  70 
c>=PP+1 

I F  ( F‘P  .  LE.  3.  )  GO  TO  20 
PP=1. 

30  X  DEL=TSAM/20.  /PP 
IF(  XDEL  .  LT.  DELT)  GO  TO  60 
PP=PP+1. 

-0  TO  30 

50  X  DEL=TSAM/ 1 0. 

60  NSAMPL=NSAMPL*DELT/XDEL 
DELT=XDEL 
GO  TO  100 

70  DELT =DELT/PP 
N'SAMPL=NSAMPL*PP 
100  CONTINUE 

* 

*  HERE  FOR  NORMAL  OPERATION 

* 

I F  ( ABS  ( KT *DELT  -  NCNT*TSAM)  .  GT.  1.  E-6  )  GO  TO  300 
IF (IT  .  EQ.  ITMAX)  NCNT=NCNT+1 
U< 1 )  =  U1 
YX  =  0. 

DO  110  1=1, NTERM 

110  YX*YX  +  ZK*ZA< I )*U( I ) 

DO  120  1=2, NTERM 
120  YX=YX  -  ZB( I )*Y( I ) 

Y< 1)  =  YX/ZB< 1 ) 

IF (  IT  .  LT.  ITMAX)  GO  TO  300 
DO  130  I=NTERM, 2, -1 
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U<I>  =  U(I-l) 

130  Y( I )  =  Y(I-l) 

* 

300  DIGTAL  ■  Y(l) 

RETURN 

END 

SUBROUTINE  RANGE 
CHARACTER  A (25) » HEAD*80 
COMMON  OUT,  XVZ (6,  402) ,  ICNT 

COMMON  /0UTP2/  1 1,  JJ,  A,  IS,  AMX (25) ,  AMN(25) ,  HEAD.  NPLOT ,  TKPLOT,  PPLANE 
COMMON  /0UTP1/  V<25),  IP,  KASE,  KS,  TSAV,  ICN 
* 

IF  (ICNT  .  GT.  1)  GO  TO  20 

DO  10  1  =  1,25 

AMX  ( I )  =  -1.  E+6 

10  AMN(  I )  =  1.  E6 

20  GO  TO  (30,  80),  IS 

30  DO  50  1=1, II 

AMX ( I )  =  AMAX 1 ( AMX ( I ) ,  V  ( I )  ) 

50  AMN(I)  =  AM INI ( AMN( I ) , V ( I )  ) 

70  RETURN 
* 

80  WRITE (OUT, 900) 

DO  100  1=1,11 

100  WRITE ( OUT ,  910)1,  A( I ) ,  AMX ( I ) , AMN ( I ) 

GO  TO  70 

900  FORMAT ( 1H0, "  MAX  AND  MIN  VALUES") 

910  FORMAT ( 1H  ,  12,  3X,  AS,  "  MAX  =",G12.  5,  "  MIN  =",G12.  5) 

END 

REAL  FUNCTION  LIMIT(U1, X,  U2) 

LIMIT  =  X 

IF  (  X  .  GT.  U2 )  LIMIT=U2 
IF  (  X  .  LT.  U1 )  LIMIT=U1 
RETURN 
END 

SUBROUTINE  TEKPLT 
* 

*  THIS  SUBROUTINE  PREPARS  A  FILE  FOR  PLOTTING 

*  ON  THE  TEKTRONICS  GRAPHIC  TERMINALS. 

*  YOU  MUST  OPEN  A  FILE  IN  YOUR  CATALOG 

* 

# 

CHARACTER  A ( 25 ) , HEAD*80 

COMMON  OUT,  X  (6,  402),  ICNT,  NVAR,  T,  TO,  TMAX,  DELT,  HDELT,  IC,  NIC,  IT,  NIT,  ITMAX,  KE 
MAI  100),  Y1PZ(  100),  Y4PZ(  100),  X1PZI  100),  X4PZI  100),  N,  KP,  NSAMPL,  KT 
COMMON  /0UTP2/  II, JJ, A,  IS, AX ( 25) , AZ ( 25 ) ,  HEAD,  NPLOT,  TKPLOT,  PPLANE 
# 

« 

IF (PPLANE  .  NE.  0)G0  TO  70 

DO  50  1=1, NVAR 

K=0 
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T=TO 

DO  20  J=l,  ICNT 
WRI TE<3, 500 ) T , X( I, J) 

K=K  ♦  1 
T  =  K*NSAMFL 
T  =  T*DELT  +  TO 
20  CONTINUE 
WRITE <3, 510) 

50  CONTINUE 
RETURN 

70  DO  80  J=l,  ICNT 

30  WRITE <3, 500 ) X < 1  * J> # X<2,  J> 

WRITE <3, 510) 

RETURN 

500  FORMAT  (El  2.  4,  El  2.  4) 

510  FORMAT  <  "  1.  E37  ,  1.  E37  " ) 

END 

REAL  FUNCTION  PTRAIN(Tl) 

COMMON  ZOUT, XVZ<6, 402),  ICNT, NVAR, T, TO,  TMAX,  DELT 
SAVE  P,  J 
DATA  J/0/ 

PTRAIN=0. 

P=.J 
F =P*T1 

IF  <  T  .  GT.  P)J=J+1 

IF< ABS(T-P)  LT.  0.  5*DELT)PTRAIN=1/BELT 

RETURN 

END 


! 


APPENDIX  I 

ANALOG  CIRCUIT  DIAGRAMS 

• 

This  appendix  contains  the  necessary  analog  flow  diagrams 
to  implement  the  Harrier  airframe  equations,  present  SAS, 
continuous  state  feedback,  and  digital  state  feedback  controller. 
Component  numbers  reference  the  present  implementation  of  the 
simulation  on  the  EAI-681  large  scale  analog  computer  in  the 
hybrid  technology  lab  of  Nimitz  Hall,  U.  S.  Naval  Academy. 
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APPENDIX  J 

DIGITAL  HYBRID  INITIATION  ROUTINE 
("AHEAD") 

Also  named  AHEAD,  this  program  sets  the  pots  and  initializes 
the  analog  computer  simulation  parameters.  Upon  completion  of 
this  task,  it  automatically  loads  the  digital  state  feedback 
control ler-SFBCON  (see  Appendix  K). 


I 

I 
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C  PROGRAMMER:  R.  V.  WALTERS 
C  TITLE:  AHEAD 

C  FUNCTION:  INITIALIZATION  OF  AV-SB 
C  HVBRID  SIMULATION  PARAMETERS 

C 

REAL  LONSTK 
REAL  LATSTK 


DIMENSION  AALON < 16,  S>,  BBLONO,  4  >,  AAL AT (11..  4>,  BBLAT  (6 
DIMENSION  AL0NC16,  1),  BLOW <12..  1>,  ALAT <16..  1>,  BLAT<1£,  1 
D I  MENS I ON  XKLON <  S ) , XKLAT < 8 > 

D I  MENS  I  ON  V5  <«,!>,  V3  <  4,  1  > 

DIMENSION  F I  LEI  <  2  >  .■  FILE2<2>..  FILE3<2>,  FIL£4<2> 
DIMENSION  FILE5<2>,  FILES<2>,  FILE?<2> 

DIMENSION  POTS < 118 > 

COMMON  XKLON..  XKLAT 

COMMON  /BLK1/VCAP, PCAP,  RCAP, PHI  CAP,  AVCAP 
COMMON  /BLK2/UCAP, WCAP, THACAP,  QCAP 
COMMON  / ELKS / D A S C A P , DRPCAP, DESCAP,  DTCAP 


DATA 

F I LEI <1 > , FI  LEI <  2 

■ /•'  A  LON 

SRC 

DATA 

FILE2<1> 

,  FILES <2 

•/•'BLON 

SRC 

DATA 

FILE3’<1> 

,  FILES <2 

■/'ALAT 

SRC 

DATA 

F I L  E  4  <  1  > 

,  F I LE4 < 2 

•X'BLAT 

SRC 

DATA 

FILES <1>,  FILES <2. 

KLAT1  , ' 

SRC 

DATA 

F I LE6  < 1 > 

,  FILES <2 

•/' KLON1 ' , ' 

SRC 

DA  TA 

F  I L  E  7  (.  1  > 

,  FI  LET’ <2 

>/' KL0N2' , ' 

SRC 

C  SET  CONSTANTS 
C 

IT-2 

CALL  H I N IT  < I E  > 

CALL  SAMO  < 1,  I  El > 

CALL  SLMCKl,  I E2  > 

CALL  STCO<0, IEZ> 

WRITE <4, 81 > 

WRITE <4, 18 > 

C 

31  FORMAT < IX,  ' ANALOG  SET  UP'  > 

82  FORMAT < IX..  •'VS  VECTOR  CALCULATED'  > 

33  FORMAT  < IX..  '  V2  VECTOR  CALCULATED  '  > 

83  FORMAT', IX,  'FILES  OPENED'  > 
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SET  SCALE  FACTORS 

VC  RP =50. 
PCAP=2.  5 
RCAF=2.  5 
PH  ICAP=2.  5 
DASCAF-IO. 
DRPCAP=10. 

A  VC  AP  =  100. 


UCAP=50. 

W  CAP =50. 

THACAP=:L.  5 
GCAP=2.  5 
DESCAP=20. 
DTCAP=20. 

INPUT  FLIGHT  VELOCITY 


WRITE (4, SI) 

FORMAT ( IX,  '  INPUT  U  (VEL  OF  A/C  C.  G. 

READ <4,  >V 
WRITE  (4,  18  > 

INPUT  PRINTING  INFO 

WRI  TE  <4,  SEC- 

FORMAT  (.  IX,  '  PR  I  NT  MA  TR I CES?  ( VES=1  >  '  > 

READ  <.4,  >  MAT  OUT 
WRITE  (.4,  18) 

WR I  TE  4,  33 

FORMAT  . IX,  PRINT  POT  SETTINGS?  <YES=1>'> 
READ  <.4,  ^  POT  OUT 
WRITE (4,  18) 

CHOOSE  LONGITUDINAL  SFB  MATRIX 
WRITE <4, 24) 

FORMATtlX,  ■'  USE  IMPROVED  KLON?  <VES=1>') 
READ <4, )KIMP 
WRITE  <4,  18) 
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SET  INITIAL  CONDITIONS  (.UNSCflLED 

VZER=0. 

F'Z£R~0. 

RZER=0. 

PH IZER-0 

UZER-0. 

WZER-0. 

THAZER-0. 

QZER--0. 

SET  INPUTS 

LGNSTK-0. 

LHTSTK— 0 
RUDFED-U. 
i  HRO 1—0. 

ALCULATE  THE  Vo  VECTOR 

VS*. 6,  i>=l. 

DO  4  I -S,  1,-1 

vs  a  ,  i>-vs<  i+i,  i>+v 

CONTINUE 
wr  i  te  4,  32  :• 


ALCULATE  THE  VS  VECTOR 

vs  <4,  i;<=i. 

DO  5  I=5:,  1  ■  -1 

VS <  I,  1  --VSU+1,  1>  +  v 

CONTINUE 

WRITE <  4,  SS> 
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C  REA  D 
C 

C 

C 


fc. 

c 

c 


AALON 


CRLL  SEEK Cl,  FI >_£1 ) 

WR I TE C 4, 88 > 

DO  6  J=l,  6 
DO  7  1=1,16 
REFID <1,  > AALON C  I,  J 
CONI I NUE 
CONTINUE 
ENDF-ILE  1 


C  EERD  BEL ON 
C 

CRLL  SEEK  C 1,  F I LE2 > 
DO  8  J=l,  T 
DO  9  1=1,8 
READ  Cl,  >  BEL ON ( I . J  > 
8  CONTINUE 

8  CONTINUE 

ENDFILE  1 
C 
C 

C  READ  RALfl T 

C 


CRLL  SEEK (1, FILES > 
DO  10  J=l,  4 
DO  11  1=1,11 
READ C 1,  )  RRLRT ( I,  J> 
11  CONTINUE 

10  CONTINUE 

ENDFILE  1 
C 
C 

C  RERD  BELAT 


CALL  SEEK Cl,  FILE4) 
DO  12  J=l,  4 
DO  IS  1=1,6 
READ  Cl,  )  BELAT C I , J > 
IS  CONTINUE 

12  CONTINUE 

ENDFILE  1 


C  READ  XKLAT 

C 

CALL  SEEK < 1,  FILES > 

DO  27  1=1..  8 
READ  Cl.-  XKLAT (  I  > 

2?  CONTINUE 

ENDFILE  1 
C 

READ  XKLON 

C 

IF  CK  IMP.  NE.  1>  GO  TO  25 
CALL  SEEKtl,  F  IL£7> 

DO  27  1=1,  3 
READ  Cl,  JXKLONC I > 

2?  CONTINUE 

ENDFILE  1 
G  Ci  T  O  2  6 

35  CONTINUE 

CALL  SEEK (1, F I LE6  > 

DO  28  I =1, 8 
READ  ( 1,  > XKLON  C I  ;> 

28  CONTINUE 

ENDFILE  1 

c 

36  CONTINUE 


C  CALCULATE  ALON 

C 

WR I  TE  c  4,  14  > 

14  FORMAT < IX,  DATA  READ' > 

C 

CALL  MMULT  ( ALON,  AALON,  16,  6,  1 
IF  CV.  LT.  5.  >  ALON <3,  l>=-32  0 
IF  CV.  LT.  5.  '?  ALON < 7,  1  >  =-3.  37 
C 

C  CALCULATE  BLON 
C 

CALL  MMULT C E-LON,  BBLON,  9,  4,  1, 
BLON C 12,  l)=BL0Nc9, 1> 

BLON C 11,  1 > =BLON C  8, 1 ) 

BLON  CIO,  1  >  =6LON  C  7,  1  > 

BLON  C  7,  1>=0. 

BLON  C  8,  1  >  =0. 

BLON  C  9,  1  >  =0 
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C  CALCULATE  ALAT 
C 

CALL  MMULTCALAT,  AALAT,  11..  4,  1,  VO 
DO  15  1=12, IE 
ALAT <  I,  1)=0. 

15  CONTINUE 
ALAT  <14,  1)=1. 

C 

C  CALCULATE  BLAT 
C 

CALL  MMULT (BLAT,  BBLAT,. 6,  4,  1,  V3) 

E;L  AT  (  S,  1  =6LAT  (.  6,  1 ) 

BLAT  1 7,  1  )  =BLAT  (.  5,  1) 

BLAT  <5,  1  )  =BLAT 4 ,  1) 

BLAT  (.4,  1  )  =BLAT  <  3,  1) 

BLAT  (3,  1  >  =0. 

BLAT  <6,  1)=0. 

DO  16  1=3,12 
BLAT  (I,  1)=0. 

16  CONTINUE 
C 

WRITE <4, 17) 

17  FORMAT  <  IX,  COEFFICIENTS  CALCULATED'  ) 

WRITE (4,  IS) 

WRITE (4, IS) 

18  FORMAT (1H  ) 


C 

13 


20 
cl'  J. 


23 


30 


IF CMATOUT.  NE.  1)  GO  TO  30 
WR I TE  <.  4,  1.9  ) 

FORMAT  C.LH1,  7X,  '  A  LON' ,  7X,  ALF.T'  ;« 
WRITE <4, IS) 

DO  20  1=1,16 

WR  I  TE  <  4,  21 )  ALON  C  1 ,  1  > ,  ALAT  <  1,1) 
CONT I NUE 

FORMAT <1X,  4X,  Ell.  4,  4X,  Ell.  4  ) 
WRITE (4, IS) 

WRITE <4, 18) 

WRITE (4, 22 > 

FORMAT  CSX,  ‘ BLOW  , ?X,  ' BLAT'  ) 
WRITE <4, IS) 

DO  23  1=1,12 

WR  I  TE  C  4,  21 )  BLON  C  1 , 1  >,  BLAT  C  1,1) 

CONTINUE 

WRITE (4, IS) 

WRITE <4, IS) 

WR  I  TE  c.  4,  23 ) 

FORMAT  <  SX,  ■'  XL. LON  ,  6X,  '  XKLAT  '  ;• 
WRITE <4, 18) 

DO  30  1=1,8 

WR  I  TE  (  4,  21 ) XKL.ON  C I  ) ,  XKL  A  T  C I  ) 

CONTINUE 

WR I  TEC 4, 18) 

WRI TEC 4, IS) 
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ZERO  THE  POTS 

DO  24  1=1,  US 
FOTS< I  )=5. 

4  CONTINUE 


CALCULATE  THE  POT  SETTINGS 

P  0  T  S  <  1  )  =  X  K  L  0  N  C  2  ;•  *  W  CAP  /  D  E  S  C  A  P 
P  0 T  S  < 2 ) = X K L  0 N  < 3 > * THACAP/DES C A P / i 0 . 
POTS  <  3  >  =XKLON «.  4  >  *GCAPXOESCAPX10. 
POTS < 4  > =XKLON < 5  >  +UCAPXDTCAP 
POTS < 5 >  =XKLON < 6 >  ♦:WCAF'XOTCAPX10. 

POT  S  <  6  >  =XKLON  <  7  )  +THACAPXDT  CAPX10. 

P  0  T  S  (  7  >  =  X  K  L  0  N  <  8  )  •+■  Q  C  A  P  /  D  T  C  A  P 
POTS < 8 ) =BLON < 5, 1 > +DT  CAP/WCAP 
P  0  T  S  <  8  )  =  X  K  L  A  T  ( 1  )  +:  V  C  A  P  /  D  A  S  C  A  P 
POTS <10 > =XKLAT <  2 ) +PCAPXDASCAPX10. 
POTS < 11 >  =L ATSTKXDASCAP 
P  0  T  S  <  1 2  >  =  X  K  L  A  T  C  3  >  +■  R  C  A  P  /  D  A  S  C  A  P 
POTS < 13 > =XKL AT <  4 ) +PH I CAPXDASCAPX10 
P  0  T  S  <  1 7  )  =  X  K  LAI  (.  5  )  *  V  C  A  P  /  D  R  P  C  A  P  X 1 8 . 

P 0 T S  < 1 8 > = X KLAT < 6 > * P C A P / D R P  C A P 
P  0  TS  <  1 8  >  =  X  K  L  A  T  <  7 )  R  C  A  P  /  D  R  P  C  A  P ,  ’10. 

P 0 TS  <2 0 )  =  BL A T <  1 ,  1  >  ♦  D A S C A P / V C A F 
POTS < 21 ) =BLAT < 2, 1 > +DRPCAP/VCAP 
P0TS<23>=.  774 
F  0  T  S 24)=.  €•  0  5 

F 0  T S  <.  2 5  >  =  6.  73 P C A P / D ASCAP/10. 

POTS  < 27  > =XKLAT < S ) +PH I CAP/DRPCAP/10 
POTS < 30 > =XKLON < 1 ) +  U C A P / D E  S C A  P  X 1 0 . 
POTS < 31 > =ALON < 5, 1 ) + UCAPXWCAP 
POTS <  32  >  =ALON (1,1) 

P  0  T  S  <  3  3  >  =  A  L  0  N  <  6 ,  1  >  X 1 0 . 

P  0  T  S  <  3  4  )  =  A  L  0  N  <  7 ,  D  +  THACA  P  X  W  C  A  p 
P 0TS<37  >=AL 0 N  < 8 ,1>*GCAPXWCAPX10. 
POTS  < 40 ) =15.  86+QCAPXDESCAPX10. 

POTS  (.42)  =  .  135 
POTS C  50 ) =ALAT C 1, 1 ) 

F'U  TS  (  51 )  =ALAT  14,  1 )  f-:  CAPXF'H  I  CAF'.-'iO 
P  0  T  S  <  5  2  >  =  A  L  A  T'  <2,1 )  *  P  C  A  P  X  V  C  A  P 
POTS  < 53 ) =ALAT <  3,  1 > + R  C  A P X V C  A P  X 1 0 . 
dotc ,  •=;.=;  ftT(4i  ,  ♦PHICAPXVCAPXIO 
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c 

c 

50 


POTS  <57)  -ALAT  <  9  ,  1  >  +VCAP/RCAP 
POTS  <  60  >  =ALON  <  2..  1  >  +WCAP/UCAP 
P  0  T  S  <  6 1 )  =  A  L  0  N  <  IS..  1  >  +  U  C  A  P  /  C  CAP 
P  0  TS<6  2  >  =  PI  L  0  N  <  4,  1  >  *  Q  C  HP/ U  C  A  P/ 1  0 . 
POTS  <  6S  )  — PlLON  < S..  1 )  +  T  H  A  C  A  P  /  IJ  C  A  P  '  1 0 . 
POTS < 64  )  =BLON <  4, 1 ) +DESC AFV WCAP 
POTS  <  67  >  =ALDN  < 12,  1  >  ♦QCAP/THACAP/iO 
P  0  T  S  <  6  9  >  -  B  L  ON  <  1 ,  1  >  *  DES  C  ft  P  /  LI  C  A  P 
PO  TS  C  70  >  -15.  8 6 + G C A P /D £ S C A P / 4 O . 

POT S  <  71  )  — .  2 4  ♦  C> A  S  C  A P / D R P CAP 
P  0  T  S  <  7  2 )  =  1 .  72 R C A P / D R P C AP/10. 

POTS <  75)=.  32 

POTS t 77  >=THROT/PTCAP 

POTS.  (.  79  )  —  E'L A  T <  4..  1 )  +DASCAP/FCAP/10. 

POTS C SO > -BLAT  < 5, 1 >  +DRPCAP/PCAP 

POTS  <.  81  >  =ALA  T  L  5,  1 )  +VCAP/PCAP 

POTS  <  83  > -ALAT  < 7,  1  >  *RCAP/PCAP/10. 

POTS  <  85  >  =ALAT  <  11,  1 )  /10. 

POTS  ( 87  > -ALAT  (.10,  1)  *  P  C  A  P  /  R  C  A  P 
P  0  T  S  ( 9  O  )  =  A  L  0  N  < 1 4 ,  1  >  f  W  C  A  P  /  0  C  A  P 


POTS <  91 )  =ALON < 16, 1 > 


P  0  T  S  <  9  7 )  =  B  L  0  N  <  1 0,  1  >  *  D  £  S  C  A  P  /  CJ  CAP  / 1 0 
POTS < 100 ) = L  0 N S  T K / D E  S C A P 
POTS <101 ) =RUDPEP/DRPCAP 


POTS  <102 ) =BLON <  2, 1 >  +DTCAP/UCAP 
P  0  T  S  <  1 0  5  )  =  B  L  0  N  <  1 1 ,  1  >  *  £>  T  C  A  P  /  Q  C  A  P  / 1 9 
P  0  T  S  <  1 1 0  >  =  A  L  A  T  <  6 ,  1  ) 


POTS <  112)=BLAT  <7,  D+DASCAP/RCAP 
POTS < 115  > =BLAT <  S, 1 ) +DRPCAP/RCAP 


WRITE <4, 50) 

FORMAT < IX,  'POT  SETTINGS 
WRITE <4, 18) 

WRITE <4, 18) 


CALCULATEO 
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C  SET  THE  POTS 

C 

CALL  SRM0(7,  I E 6 > 

C 

DO  51  1=1,  11 8 
ICOF=ABS(10O00.  *POTS(  I  >  > 

IF  (ICOF.  EQ  50000;-  GO  TO  51 
IF  (ICOF  LE.  10fi00 >  go  TO  25 
WRITE (.4,  2E >  I 

26  FORMAT (IX, ' POT  SETTING  TOO  LARGE 

I C OF -9995 
25  CONTINUE 

CALL  SPOT  (I,  ICOF,  IT,  IE7> 

C 

IF (POTOUT.  NE.  1>  GO  TO  51 

C 

WRITE <4,  >  POT'S  <  I  >,  IE? 

51  CONTINUE 


C 

C  ALATC9, 1>  SIGN  TEST 

C 

IF  CALATC9,  1 J.  GT.  O.  >  GO  TO  52 
CALL  SSRF’dl,  I E4  > 

GO  TO  110 

52  CALL  SSRM<11, IE5) 

110  CONTINUE 

C 

c 


CALL  SAM Ci ( 1,  I E 8 > 

C 

WR I TE  C 4,  18 > 

WR I  TE  c  4,  54  :• 

WRITE <4,  18 > 

54  FORMAT (IX, 'POT  SETTING 

C 

C  BEGIN  SIMULATION 

C 

c 

CALL  SFBCON ( V) 

C 

STOP 

END 


COMPLETE 
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APPENDIX  K 

DIGITAL  STATE  FEEDBACK  CONTROLLER  SUBROUTINE 
("SFBCON") 

This  program  implements  the  digital  state  feedback  controller. 
Note  that  it  operates  continuously  in  an  endless  loop.  It  is 
activated  from  the  logic  and  switching  controls  of  the  analog 
computer.  It  has  been  determined  that  its  sampling  rate  is 
approximately  71  samples  per  second  (T  =  14  ’’s.). 
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SUBROUT I NE  SF6CGN ( V > 


D I  MENS  I  ON  XKLON  C  8  ) ,  XKLAT  ( 8  > 

COMMON  XKLON, XKLAT 

COMMON  /BLK1/VCAP, PCAP, RCAP, PH  I  CAP, AVCAP 
COMMON  /BLK2/UCAP, WCAP, THACAP, QChP 
COMMON  / B L K 3 / D  A  S  C  A  P , DRPCfiP, DESCAP, DTCAP 

WR I T  E ( 4, 2> 

FORMAT (IX, ' DIGITAL  CONTROL  ACTIVATED') 
CONTINUE 

READ  AND  CONVERT  ANALOG  STATES 

CALL  CRACS <8,  IV,  I  El) 

CALL  CRACS ( 7, IP, IE2) 

CALL  CRACS <8,  IR,  IE2) 

CALL  CRACS (9, I PH I, IE4) 

CALL  CRACS <2, IU, IE5) 

CALL  CRACS <3, IW, IE6) 

CALL  CRACS (4, ITHA, IE7) 

CALL  CRACS (5, 10, IES) 


CALC  LATERAL  FEEDBACK  SIGNAL  (SCALED) 

I DAS*.  1* ( XKLAT ( 1 ) *  I V+VCAP+XKLAT ( 2 ) *  I P*PCAP ) XDASC AP 
I DAS* I DAS+.  1* C  XKLAT (3 >  *  I R+RCAP+XKLAT ( 4 ) *  I  PH  I *PH I  CAP) /DASCA 
I DAS*- IDAS 

I DRP*.  1+ (  XKLA T ( 5 ) *  I V+  VCAP+XKLAT ( 6 ) + I P+PCAP ) XDRPCAP 
I DRP* I DRP+.  1* C  XKLAT ( 7 ) +RCAP+ 1 R+XKLAT ( 8 ) *  I  PH  I *PH I  CAP )  XDRPCA 

CALC  LONGITUDINAL  CONTROL  INPUT  (SCALED) 


IDES*.  1+  C XKLON C 1)+  I U+UCAP+XKLON ( 2 ) *  I W+  WCAP ) /DESCAP 

I DES* I DES+.  1* ( XKLON ( 3 > *  I THA+THACAP+XKLON ( 4  >  *  I Q+QCAP ) XDESCA 

I DT*  <  XKLON ( 5 ) *  I U+UCAP+XKLON ( 6 ) *  I U+  WCAP ) /DTCAP 

IDT*  I DT+  ( XKLON  (  7 ) *  I THA+ THACAP+XKLON ( 8 )  * 1 GJ+QCAP >/DTCAP 

TEST  FOR  SATURATION 

IF( IDAS.  GT  10OO0)  IDAS*10000 
IF  ( I  DRP.  GT.  10000)  I  DRP*10000 
IF ( IDES  GT  10000)  I DES*10000 
I F <  IDT,  GT  10000)  IDT *10000 


c 

c 


c 

c 

c 


PERFORM  DIGITAL  TO  ANALOG  CONVERSION 

CALL  LTDASC0, IDAS, IE9) 

CALL  TLDA 

CALL  LTDASC1, IDRP, IE10) 

CALL  TLDA 

CALL  LTDASC2, IDES, IEll) 

CALL  TLDA 

CALL  LTDASO,  IDT,  IE12) 

CALL  TLDA 


GO  TO  1 

RETURN 

END 


SUBROUTINE  MMULTCR,  A,  N,  M,  L,  B> 

D I  MENS  I  ON  R  C  N,  L  > ,  A  C  N,  M ) ,  B  C  M,  L  > 
DIMENSION  VC 20, 20) 

C 

DO  101  1=1,  N 
DO  iOl  J=l,  L 
VC  I,  J>=0. 

DO  101  K=l,  M 

V  < I ,  J  >  =A  C I ,  K > +B <  K, J >  +V < I , J  > 

101  CONTINUE 

DO  102  1=1,  N 
DO  102  J=i,  L 
RC  I,  J)=VC  I,  J) 

102  CONTINUE 
RETURN 
END 
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APPENDIX  L 

AV-3B  OUTPUT  TIME  TRACES 

This  appendix  contains  time  traces  of  the  AHEAD  AV-8B, 
illustrating  the  results  of  the  research.  In  all  cases,  the  x-axis 
variable  is  time  in  seconds,  and  the  ordinate  represents  the  system 
states  in: 

Radians  -  0,  <J> 

Radians  per  sec  -  Q,  P,  R 
Feet  per  sec  -  U,  W,  V 

The  results  are  divided  into  six  sections: 

I.  V  =  30;  No  stability  augmentation;  no  I.C.'s. 

Shows  open  loop  system  response  to  specified  control 

inputs.  Simulation  fidelity  can  be  assessed  by  comparison 

with  Cal  span  results,  Aopendix  M. 

II-IV 

<J>0  =  15°  =  .26R 
vQ  =  10  kt  -  16.9  fps 
rQ  =  5°/sec  -  .09  R/s 
eQ  =  15°  -  .  26R 

Initial  condition  response  for  sections  II  through  IV 
is  meant  to  be  representative  of  the  onset  of  the  coupled 
intake  momentum  drag  pehnomena. 

II.  LATERAL  RESPONSE;  V  =  0,  30,  50,  80 

Compares  performance  of  present  SAS  and  state  feedback 
controller  at  various  velocities. 
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APPENDIX  L 


f, 

44 

(Continued) 

III.  LONGITUDINAL  RESPONSE;  V  =  0,  30,  50,  80 
Same  as  II,  except  longitudinal  case. 

IV.  LATERAL  AND  LONGITUDINAL;  V  =  30 

Compares  performance  of  the  various  stability 
augmentation  options. 

V.  LONGITUDINAL;  V  =  30;  INPUT  PITCH  ANGLE 

Shows  response  of  SFB-controlled  system  to  0Q  = 

15°,  20°,  25°,  and  30°. 

VI.  LATERAL;  rQ  =  5°/s;  4»0  =  15°;  INPUT  SIDE  VELOCITY 

Shows  response  of  SFB-controlled  system  to  Vo  =  10  » 
20,  30  and  40  kts.  in  presence  of  an  initial  bank 

angle  and  yaw  rate.  (V*3o) 
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APPENDIX  M 

CALSPAN  TIME  RESPONSES 
(EXTRACTED  FROM  TM  NO.  98) 

This  appendix  contains  Calspan  results  to  control  inputs 
used  in  Appendix  L,  Part  I.  Comparison  of  the  two  enables  assess¬ 
ment  of  simulation  fidelity  to  be  made. 
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